| // Written in the D programming language. |
| |
| /** |
| This is a submodule of $(MREF std, math). |
| |
| It contains several functions for work with floating point numbers. |
| |
| Copyright: Copyright The D Language Foundation 2000 - 2011. |
| License: $(HTTP www.boost.org/LICENSE_1_0.txt, Boost License 1.0). |
| Authors: $(HTTP digitalmars.com, Walter Bright), Don Clugston, |
| Conversion of CEPHES math library to D by Iain Buclaw and David Nadlinger |
| Source: $(PHOBOSSRC std/math/operations.d) |
| |
| Macros: |
| TABLE_SV = <table border="1" cellpadding="4" cellspacing="0"> |
| <caption>Special Values</caption> |
| $0</table> |
| SVH = $(TR $(TH $1) $(TH $2)) |
| SV = $(TR $(TD $1) $(TD $2)) |
| NAN = $(RED NAN) |
| PLUSMN = ± |
| INFIN = ∞ |
| LT = < |
| GT = > |
| */ |
| |
| module std.math.operations; |
| |
| import std.traits : CommonType, isFloatingPoint, isIntegral, Unqual; |
| |
| // Functions for NaN payloads |
| /* |
| * A 'payload' can be stored in the significand of a $(NAN). One bit is required |
| * to distinguish between a quiet and a signalling $(NAN). This leaves 22 bits |
| * of payload for a float; 51 bits for a double; 62 bits for an 80-bit real; |
| * and 111 bits for a 128-bit quad. |
| */ |
| /** |
| * Create a quiet $(NAN), storing an integer inside the payload. |
| * |
| * For floats, the largest possible payload is 0x3F_FFFF. |
| * For doubles, it is 0x3_FFFF_FFFF_FFFF. |
| * For 80-bit or 128-bit reals, it is 0x3FFF_FFFF_FFFF_FFFF. |
| */ |
| real NaN(ulong payload) @trusted pure nothrow @nogc |
| { |
| import std.math.traits : floatTraits, RealFormat; |
| |
| alias F = floatTraits!(real); |
| static if (F.realFormat == RealFormat.ieeeExtended || |
| F.realFormat == RealFormat.ieeeExtended53) |
| { |
| // real80 (in x86 real format, the implied bit is actually |
| // not implied but a real bit which is stored in the real) |
| ulong v = 3; // implied bit = 1, quiet bit = 1 |
| } |
| else |
| { |
| ulong v = 1; // no implied bit. quiet bit = 1 |
| } |
| if (__ctfe) |
| { |
| v = 1; // We use a double in CTFE. |
| assert(payload >>> 51 == 0, |
| "Cannot set more than 51 bits of NaN payload in CTFE."); |
| } |
| |
| |
| ulong a = payload; |
| |
| // 22 Float bits |
| ulong w = a & 0x3F_FFFF; |
| a -= w; |
| |
| v <<=22; |
| v |= w; |
| a >>=22; |
| |
| // 29 Double bits |
| v <<=29; |
| w = a & 0xFFF_FFFF; |
| v |= w; |
| a -= w; |
| a >>=29; |
| |
| if (__ctfe) |
| { |
| v |= 0x7FF0_0000_0000_0000; |
| return *cast(double*) &v; |
| } |
| else static if (F.realFormat == RealFormat.ieeeDouble) |
| { |
| v |= 0x7FF0_0000_0000_0000; |
| real x; |
| * cast(ulong *)(&x) = v; |
| return x; |
| } |
| else |
| { |
| v <<=11; |
| a &= 0x7FF; |
| v |= a; |
| real x = real.nan; |
| |
| // Extended real bits |
| static if (F.realFormat == RealFormat.ieeeQuadruple) |
| { |
| v <<= 1; // there's no implicit bit |
| |
| version (LittleEndian) |
| { |
| *cast(ulong*)(6+cast(ubyte*)(&x)) = v; |
| } |
| else |
| { |
| *cast(ulong*)(2+cast(ubyte*)(&x)) = v; |
| } |
| } |
| else |
| { |
| *cast(ulong *)(&x) = v; |
| } |
| return x; |
| } |
| } |
| |
| /// |
| @safe @nogc pure nothrow unittest |
| { |
| import std.math.traits : isNaN; |
| |
| real a = NaN(1_000_000); |
| assert(isNaN(a)); |
| assert(getNaNPayload(a) == 1_000_000); |
| } |
| |
| @system pure nothrow @nogc unittest // not @safe because taking address of local. |
| { |
| import std.math.traits : floatTraits, RealFormat; |
| |
| static if (floatTraits!(real).realFormat == RealFormat.ieeeDouble) |
| { |
| auto x = NaN(1); |
| auto xl = *cast(ulong*)&x; |
| assert(xl & 0x8_0000_0000_0000UL); //non-signaling bit, bit 52 |
| assert((xl & 0x7FF0_0000_0000_0000UL) == 0x7FF0_0000_0000_0000UL); //all exp bits set |
| } |
| } |
| |
| /** |
| * Extract an integral payload from a $(NAN). |
| * |
| * Returns: |
| * the integer payload as a ulong. |
| * |
| * For floats, the largest possible payload is 0x3F_FFFF. |
| * For doubles, it is 0x3_FFFF_FFFF_FFFF. |
| * For 80-bit or 128-bit reals, it is 0x3FFF_FFFF_FFFF_FFFF. |
| */ |
| ulong getNaNPayload(real x) @trusted pure nothrow @nogc |
| { |
| import std.math.traits : floatTraits, RealFormat; |
| |
| // assert(isNaN(x)); |
| alias F = floatTraits!(real); |
| ulong m = void; |
| if (__ctfe) |
| { |
| double y = x; |
| m = *cast(ulong*) &y; |
| // Make it look like an 80-bit significand. |
| // Skip exponent, and quiet bit |
| m &= 0x0007_FFFF_FFFF_FFFF; |
| m <<= 11; |
| } |
| else static if (F.realFormat == RealFormat.ieeeDouble) |
| { |
| m = *cast(ulong*)(&x); |
| // Make it look like an 80-bit significand. |
| // Skip exponent, and quiet bit |
| m &= 0x0007_FFFF_FFFF_FFFF; |
| m <<= 11; |
| } |
| else static if (F.realFormat == RealFormat.ieeeQuadruple) |
| { |
| version (LittleEndian) |
| { |
| m = *cast(ulong*)(6+cast(ubyte*)(&x)); |
| } |
| else |
| { |
| m = *cast(ulong*)(2+cast(ubyte*)(&x)); |
| } |
| |
| m >>= 1; // there's no implicit bit |
| } |
| else |
| { |
| m = *cast(ulong*)(&x); |
| } |
| |
| // ignore implicit bit and quiet bit |
| |
| const ulong f = m & 0x3FFF_FF00_0000_0000L; |
| |
| ulong w = f >>> 40; |
| w |= (m & 0x00FF_FFFF_F800L) << (22 - 11); |
| w |= (m & 0x7FF) << 51; |
| return w; |
| } |
| |
| /// |
| @safe @nogc pure nothrow unittest |
| { |
| import std.math.traits : isNaN; |
| |
| real a = NaN(1_000_000); |
| assert(isNaN(a)); |
| assert(getNaNPayload(a) == 1_000_000); |
| } |
| |
| @safe @nogc pure nothrow unittest |
| { |
| import std.math.traits : isIdentical, isNaN; |
| |
| enum real a = NaN(1_000_000); |
| static assert(isNaN(a)); |
| static assert(getNaNPayload(a) == 1_000_000); |
| real b = NaN(1_000_000); |
| assert(isIdentical(b, a)); |
| // The CTFE version of getNaNPayload relies on it being impossible |
| // for a CTFE-constructed NaN to have more than 51 bits of payload. |
| enum nanNaN = NaN(getNaNPayload(real.nan)); |
| assert(isIdentical(real.nan, nanNaN)); |
| static if (real.init != real.init) |
| { |
| enum initNaN = NaN(getNaNPayload(real.init)); |
| assert(isIdentical(real.init, initNaN)); |
| } |
| } |
| |
| debug(UnitTest) |
| { |
| @safe pure nothrow @nogc unittest |
| { |
| real nan4 = NaN(0x789_ABCD_EF12_3456); |
| static if (floatTraits!(real).realFormat == RealFormat.ieeeExtended |
| || floatTraits!(real).realFormat == RealFormat.ieeeQuadruple) |
| { |
| assert(getNaNPayload(nan4) == 0x789_ABCD_EF12_3456); |
| } |
| else |
| { |
| assert(getNaNPayload(nan4) == 0x1_ABCD_EF12_3456); |
| } |
| double nan5 = nan4; |
| assert(getNaNPayload(nan5) == 0x1_ABCD_EF12_3456); |
| float nan6 = nan4; |
| assert(getNaNPayload(nan6) == 0x12_3456); |
| nan4 = NaN(0xFABCD); |
| assert(getNaNPayload(nan4) == 0xFABCD); |
| nan6 = nan4; |
| assert(getNaNPayload(nan6) == 0xFABCD); |
| nan5 = NaN(0x100_0000_0000_3456); |
| assert(getNaNPayload(nan5) == 0x0000_0000_3456); |
| } |
| } |
| |
| /** |
| * Calculate the next largest floating point value after x. |
| * |
| * Return the least number greater than x that is representable as a real; |
| * thus, it gives the next point on the IEEE number line. |
| * |
| * $(TABLE_SV |
| * $(SVH x, nextUp(x) ) |
| * $(SV -$(INFIN), -real.max ) |
| * $(SV $(PLUSMN)0.0, real.min_normal*real.epsilon ) |
| * $(SV real.max, $(INFIN) ) |
| * $(SV $(INFIN), $(INFIN) ) |
| * $(SV $(NAN), $(NAN) ) |
| * ) |
| */ |
| real nextUp(real x) @trusted pure nothrow @nogc |
| { |
| import std.math.traits : floatTraits, RealFormat, MANTISSA_MSB, MANTISSA_LSB; |
| |
| alias F = floatTraits!(real); |
| static if (F.realFormat != RealFormat.ieeeDouble) |
| { |
| if (__ctfe) |
| { |
| if (x == -real.infinity) |
| return -real.max; |
| if (!(x < real.infinity)) // Infinity or NaN. |
| return x; |
| real delta; |
| // Start with a decent estimate of delta. |
| if (x <= 0x1.ffffffffffffep+1023 && x >= -double.max) |
| { |
| const double d = cast(double) x; |
| delta = (cast(real) nextUp(d) - cast(real) d) * 0x1p-11L; |
| while (x + (delta * 0x1p-100L) > x) |
| delta *= 0x1p-100L; |
| } |
| else |
| { |
| delta = 0x1p960L; |
| while (!(x + delta > x) && delta < real.max * 0x1p-100L) |
| delta *= 0x1p100L; |
| } |
| if (x + delta > x) |
| { |
| while (x + (delta / 2) > x) |
| delta /= 2; |
| } |
| else |
| { |
| do { delta += delta; } while (!(x + delta > x)); |
| } |
| if (x < 0 && x + delta == 0) |
| return -0.0L; |
| return x + delta; |
| } |
| } |
| static if (F.realFormat == RealFormat.ieeeDouble) |
| { |
| return nextUp(cast(double) x); |
| } |
| else static if (F.realFormat == RealFormat.ieeeQuadruple) |
| { |
| ushort e = F.EXPMASK & (cast(ushort *)&x)[F.EXPPOS_SHORT]; |
| if (e == F.EXPMASK) |
| { |
| // NaN or Infinity |
| if (x == -real.infinity) return -real.max; |
| return x; // +Inf and NaN are unchanged. |
| } |
| |
| auto ps = cast(ulong *)&x; |
| if (ps[MANTISSA_MSB] & 0x8000_0000_0000_0000) |
| { |
| // Negative number |
| if (ps[MANTISSA_LSB] == 0 && ps[MANTISSA_MSB] == 0x8000_0000_0000_0000) |
| { |
| // it was negative zero, change to smallest subnormal |
| ps[MANTISSA_LSB] = 1; |
| ps[MANTISSA_MSB] = 0; |
| return x; |
| } |
| if (ps[MANTISSA_LSB] == 0) --ps[MANTISSA_MSB]; |
| --ps[MANTISSA_LSB]; |
| } |
| else |
| { |
| // Positive number |
| ++ps[MANTISSA_LSB]; |
| if (ps[MANTISSA_LSB] == 0) ++ps[MANTISSA_MSB]; |
| } |
| return x; |
| } |
| else static if (F.realFormat == RealFormat.ieeeExtended || |
| F.realFormat == RealFormat.ieeeExtended53) |
| { |
| // For 80-bit reals, the "implied bit" is a nuisance... |
| ushort *pe = cast(ushort *)&x; |
| ulong *ps = cast(ulong *)&x; |
| // EPSILON is 1 for 64-bit, and 2048 for 53-bit precision reals. |
| enum ulong EPSILON = 2UL ^^ (64 - real.mant_dig); |
| |
| if ((pe[F.EXPPOS_SHORT] & F.EXPMASK) == F.EXPMASK) |
| { |
| // First, deal with NANs and infinity |
| if (x == -real.infinity) return -real.max; |
| return x; // +Inf and NaN are unchanged. |
| } |
| if (pe[F.EXPPOS_SHORT] & 0x8000) |
| { |
| // Negative number -- need to decrease the significand |
| *ps -= EPSILON; |
| // Need to mask with 0x7FFF... so subnormals are treated correctly. |
| if ((*ps & 0x7FFF_FFFF_FFFF_FFFF) == 0x7FFF_FFFF_FFFF_FFFF) |
| { |
| if (pe[F.EXPPOS_SHORT] == 0x8000) // it was negative zero |
| { |
| *ps = 1; |
| pe[F.EXPPOS_SHORT] = 0; // smallest subnormal. |
| return x; |
| } |
| |
| --pe[F.EXPPOS_SHORT]; |
| |
| if (pe[F.EXPPOS_SHORT] == 0x8000) |
| return x; // it's become a subnormal, implied bit stays low. |
| |
| *ps = 0xFFFF_FFFF_FFFF_FFFF; // set the implied bit |
| return x; |
| } |
| return x; |
| } |
| else |
| { |
| // Positive number -- need to increase the significand. |
| // Works automatically for positive zero. |
| *ps += EPSILON; |
| if ((*ps & 0x7FFF_FFFF_FFFF_FFFF) == 0) |
| { |
| // change in exponent |
| ++pe[F.EXPPOS_SHORT]; |
| *ps = 0x8000_0000_0000_0000; // set the high bit |
| } |
| } |
| return x; |
| } |
| else // static if (F.realFormat == RealFormat.ibmExtended) |
| { |
| assert(0, "nextUp not implemented"); |
| } |
| } |
| |
| /** ditto */ |
| double nextUp(double x) @trusted pure nothrow @nogc |
| { |
| ulong s = *cast(ulong *)&x; |
| |
| if ((s & 0x7FF0_0000_0000_0000) == 0x7FF0_0000_0000_0000) |
| { |
| // First, deal with NANs and infinity |
| if (x == -x.infinity) return -x.max; |
| return x; // +INF and NAN are unchanged. |
| } |
| if (s & 0x8000_0000_0000_0000) // Negative number |
| { |
| if (s == 0x8000_0000_0000_0000) // it was negative zero |
| { |
| s = 0x0000_0000_0000_0001; // change to smallest subnormal |
| return *cast(double*) &s; |
| } |
| --s; |
| } |
| else |
| { // Positive number |
| ++s; |
| } |
| return *cast(double*) &s; |
| } |
| |
| /** ditto */ |
| float nextUp(float x) @trusted pure nothrow @nogc |
| { |
| uint s = *cast(uint *)&x; |
| |
| if ((s & 0x7F80_0000) == 0x7F80_0000) |
| { |
| // First, deal with NANs and infinity |
| if (x == -x.infinity) return -x.max; |
| |
| return x; // +INF and NAN are unchanged. |
| } |
| if (s & 0x8000_0000) // Negative number |
| { |
| if (s == 0x8000_0000) // it was negative zero |
| { |
| s = 0x0000_0001; // change to smallest subnormal |
| return *cast(float*) &s; |
| } |
| |
| --s; |
| } |
| else |
| { |
| // Positive number |
| ++s; |
| } |
| return *cast(float*) &s; |
| } |
| |
| /// |
| @safe @nogc pure nothrow unittest |
| { |
| assert(nextUp(1.0 - 1.0e-6).feqrel(0.999999) > 16); |
| assert(nextUp(1.0 - real.epsilon).feqrel(1.0) > 16); |
| } |
| |
| /** |
| * Calculate the next smallest floating point value before x. |
| * |
| * Return the greatest number less than x that is representable as a real; |
| * thus, it gives the previous point on the IEEE number line. |
| * |
| * $(TABLE_SV |
| * $(SVH x, nextDown(x) ) |
| * $(SV $(INFIN), real.max ) |
| * $(SV $(PLUSMN)0.0, -real.min_normal*real.epsilon ) |
| * $(SV -real.max, -$(INFIN) ) |
| * $(SV -$(INFIN), -$(INFIN) ) |
| * $(SV $(NAN), $(NAN) ) |
| * ) |
| */ |
| real nextDown(real x) @safe pure nothrow @nogc |
| { |
| return -nextUp(-x); |
| } |
| |
| /** ditto */ |
| double nextDown(double x) @safe pure nothrow @nogc |
| { |
| return -nextUp(-x); |
| } |
| |
| /** ditto */ |
| float nextDown(float x) @safe pure nothrow @nogc |
| { |
| return -nextUp(-x); |
| } |
| |
| /// |
| @safe pure nothrow @nogc unittest |
| { |
| assert( nextDown(1.0 + real.epsilon) == 1.0); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| import std.math.traits : floatTraits, RealFormat, isIdentical; |
| |
| static if (floatTraits!(real).realFormat == RealFormat.ieeeExtended || |
| floatTraits!(real).realFormat == RealFormat.ieeeDouble || |
| floatTraits!(real).realFormat == RealFormat.ieeeExtended53 || |
| floatTraits!(real).realFormat == RealFormat.ieeeQuadruple) |
| { |
| // Tests for reals |
| assert(isIdentical(nextUp(NaN(0xABC)), NaN(0xABC))); |
| //static assert(isIdentical(nextUp(NaN(0xABC)), NaN(0xABC))); |
| // negative numbers |
| assert( nextUp(-real.infinity) == -real.max ); |
| assert( nextUp(-1.0L-real.epsilon) == -1.0 ); |
| assert( nextUp(-2.0L) == -2.0 + real.epsilon); |
| static assert( nextUp(-real.infinity) == -real.max ); |
| static assert( nextUp(-1.0L-real.epsilon) == -1.0 ); |
| static assert( nextUp(-2.0L) == -2.0 + real.epsilon); |
| // subnormals and zero |
| assert( nextUp(-real.min_normal) == -real.min_normal*(1-real.epsilon) ); |
| assert( nextUp(-real.min_normal*(1-real.epsilon)) == -real.min_normal*(1-2*real.epsilon) ); |
| assert( isIdentical(-0.0L, nextUp(-real.min_normal*real.epsilon)) ); |
| assert( nextUp(-0.0L) == real.min_normal*real.epsilon ); |
| assert( nextUp(0.0L) == real.min_normal*real.epsilon ); |
| assert( nextUp(real.min_normal*(1-real.epsilon)) == real.min_normal ); |
| assert( nextUp(real.min_normal) == real.min_normal*(1+real.epsilon) ); |
| static assert( nextUp(-real.min_normal) == -real.min_normal*(1-real.epsilon) ); |
| static assert( nextUp(-real.min_normal*(1-real.epsilon)) == -real.min_normal*(1-2*real.epsilon) ); |
| static assert( -0.0L is nextUp(-real.min_normal*real.epsilon) ); |
| static assert( nextUp(-0.0L) == real.min_normal*real.epsilon ); |
| static assert( nextUp(0.0L) == real.min_normal*real.epsilon ); |
| static assert( nextUp(real.min_normal*(1-real.epsilon)) == real.min_normal ); |
| static assert( nextUp(real.min_normal) == real.min_normal*(1+real.epsilon) ); |
| // positive numbers |
| assert( nextUp(1.0L) == 1.0 + real.epsilon ); |
| assert( nextUp(2.0L-real.epsilon) == 2.0 ); |
| assert( nextUp(real.max) == real.infinity ); |
| assert( nextUp(real.infinity)==real.infinity ); |
| static assert( nextUp(1.0L) == 1.0 + real.epsilon ); |
| static assert( nextUp(2.0L-real.epsilon) == 2.0 ); |
| static assert( nextUp(real.max) == real.infinity ); |
| static assert( nextUp(real.infinity)==real.infinity ); |
| // ctfe near double.max boundary |
| static assert(nextUp(nextDown(cast(real) double.max)) == cast(real) double.max); |
| } |
| |
| double n = NaN(0xABC); |
| assert(isIdentical(nextUp(n), n)); |
| // negative numbers |
| assert( nextUp(-double.infinity) == -double.max ); |
| assert( nextUp(-1-double.epsilon) == -1.0 ); |
| assert( nextUp(-2.0) == -2.0 + double.epsilon); |
| // subnormals and zero |
| |
| assert( nextUp(-double.min_normal) == -double.min_normal*(1-double.epsilon) ); |
| assert( nextUp(-double.min_normal*(1-double.epsilon)) == -double.min_normal*(1-2*double.epsilon) ); |
| assert( isIdentical(-0.0, nextUp(-double.min_normal*double.epsilon)) ); |
| assert( nextUp(0.0) == double.min_normal*double.epsilon ); |
| assert( nextUp(-0.0) == double.min_normal*double.epsilon ); |
| assert( nextUp(double.min_normal*(1-double.epsilon)) == double.min_normal ); |
| assert( nextUp(double.min_normal) == double.min_normal*(1+double.epsilon) ); |
| // positive numbers |
| assert( nextUp(1.0) == 1.0 + double.epsilon ); |
| assert( nextUp(2.0-double.epsilon) == 2.0 ); |
| assert( nextUp(double.max) == double.infinity ); |
| |
| float fn = NaN(0xABC); |
| assert(isIdentical(nextUp(fn), fn)); |
| float f = -float.min_normal*(1-float.epsilon); |
| float f1 = -float.min_normal; |
| assert( nextUp(f1) == f); |
| f = 1.0f+float.epsilon; |
| f1 = 1.0f; |
| assert( nextUp(f1) == f ); |
| f1 = -0.0f; |
| assert( nextUp(f1) == float.min_normal*float.epsilon); |
| assert( nextUp(float.infinity)==float.infinity ); |
| |
| assert(nextDown(1.0L+real.epsilon)==1.0); |
| assert(nextDown(1.0+double.epsilon)==1.0); |
| f = 1.0f+float.epsilon; |
| assert(nextDown(f)==1.0); |
| assert(nextafter(1.0+real.epsilon, -real.infinity)==1.0); |
| |
| // CTFE |
| |
| enum double ctfe_n = NaN(0xABC); |
| //static assert(isIdentical(nextUp(ctfe_n), ctfe_n)); // FIXME: https://issues.dlang.org/show_bug.cgi?id=20197 |
| static assert(nextUp(double.nan) is double.nan); |
| // negative numbers |
| static assert( nextUp(-double.infinity) == -double.max ); |
| static assert( nextUp(-1-double.epsilon) == -1.0 ); |
| static assert( nextUp(-2.0) == -2.0 + double.epsilon); |
| // subnormals and zero |
| |
| static assert( nextUp(-double.min_normal) == -double.min_normal*(1-double.epsilon) ); |
| static assert( nextUp(-double.min_normal*(1-double.epsilon)) == -double.min_normal*(1-2*double.epsilon) ); |
| static assert( -0.0 is nextUp(-double.min_normal*double.epsilon) ); |
| static assert( nextUp(0.0) == double.min_normal*double.epsilon ); |
| static assert( nextUp(-0.0) == double.min_normal*double.epsilon ); |
| static assert( nextUp(double.min_normal*(1-double.epsilon)) == double.min_normal ); |
| static assert( nextUp(double.min_normal) == double.min_normal*(1+double.epsilon) ); |
| // positive numbers |
| static assert( nextUp(1.0) == 1.0 + double.epsilon ); |
| static assert( nextUp(2.0-double.epsilon) == 2.0 ); |
| static assert( nextUp(double.max) == double.infinity ); |
| |
| enum float ctfe_fn = NaN(0xABC); |
| //static assert(isIdentical(nextUp(ctfe_fn), ctfe_fn)); // FIXME: https://issues.dlang.org/show_bug.cgi?id=20197 |
| static assert(nextUp(float.nan) is float.nan); |
| static assert(nextUp(-float.min_normal) == -float.min_normal*(1-float.epsilon)); |
| static assert(nextUp(1.0f) == 1.0f+float.epsilon); |
| static assert(nextUp(-0.0f) == float.min_normal*float.epsilon); |
| static assert(nextUp(float.infinity)==float.infinity); |
| static assert(nextDown(1.0L+real.epsilon)==1.0); |
| static assert(nextDown(1.0+double.epsilon)==1.0); |
| static assert(nextDown(1.0f+float.epsilon)==1.0); |
| static assert(nextafter(1.0+real.epsilon, -real.infinity)==1.0); |
| } |
| |
| |
| |
| /****************************************** |
| * Calculates the next representable value after x in the direction of y. |
| * |
| * If y > x, the result will be the next largest floating-point value; |
| * if y < x, the result will be the next smallest value. |
| * If x == y, the result is y. |
| * If x or y is a NaN, the result is a NaN. |
| * |
| * Remarks: |
| * This function is not generally very useful; it's almost always better to use |
| * the faster functions nextUp() or nextDown() instead. |
| * |
| * The FE_INEXACT and FE_OVERFLOW exceptions will be raised if x is finite and |
| * the function result is infinite. The FE_INEXACT and FE_UNDERFLOW |
| * exceptions will be raised if the function value is subnormal, and x is |
| * not equal to y. |
| */ |
| T nextafter(T)(const T x, const T y) @safe pure nothrow @nogc |
| { |
| import std.math.traits : isNaN; |
| |
| if (x == y || isNaN(y)) |
| { |
| return y; |
| } |
| |
| if (isNaN(x)) |
| { |
| return x; |
| } |
| |
| return ((y>x) ? nextUp(x) : nextDown(x)); |
| } |
| |
| /// |
| @safe pure nothrow @nogc unittest |
| { |
| import std.math.traits : isNaN; |
| |
| float a = 1; |
| assert(is(typeof(nextafter(a, a)) == float)); |
| assert(nextafter(a, a.infinity) > a); |
| assert(isNaN(nextafter(a, a.nan))); |
| assert(isNaN(nextafter(a.nan, a))); |
| |
| double b = 2; |
| assert(is(typeof(nextafter(b, b)) == double)); |
| assert(nextafter(b, b.infinity) > b); |
| assert(isNaN(nextafter(b, b.nan))); |
| assert(isNaN(nextafter(b.nan, b))); |
| |
| real c = 3; |
| assert(is(typeof(nextafter(c, c)) == real)); |
| assert(nextafter(c, c.infinity) > c); |
| assert(isNaN(nextafter(c, c.nan))); |
| assert(isNaN(nextafter(c.nan, c))); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| import std.math.traits : isNaN, signbit; |
| |
| // CTFE |
| enum float a = 1; |
| static assert(is(typeof(nextafter(a, a)) == float)); |
| static assert(nextafter(a, a.infinity) > a); |
| static assert(isNaN(nextafter(a, a.nan))); |
| static assert(isNaN(nextafter(a.nan, a))); |
| |
| enum double b = 2; |
| static assert(is(typeof(nextafter(b, b)) == double)); |
| static assert(nextafter(b, b.infinity) > b); |
| static assert(isNaN(nextafter(b, b.nan))); |
| static assert(isNaN(nextafter(b.nan, b))); |
| |
| enum real c = 3; |
| static assert(is(typeof(nextafter(c, c)) == real)); |
| static assert(nextafter(c, c.infinity) > c); |
| static assert(isNaN(nextafter(c, c.nan))); |
| static assert(isNaN(nextafter(c.nan, c))); |
| |
| enum real negZero = nextafter(+0.0L, -0.0L); |
| static assert(negZero == -0.0L); |
| static assert(signbit(negZero)); |
| |
| static assert(nextafter(c, c) == c); |
| } |
| |
| //real nexttoward(real x, real y) { return core.stdc.math.nexttowardl(x, y); } |
| |
| /** |
| * Returns the positive difference between x and y. |
| * |
| * Equivalent to `fmax(x-y, 0)`. |
| * |
| * Returns: |
| * $(TABLE_SV |
| * $(TR $(TH x, y) $(TH fdim(x, y))) |
| * $(TR $(TD x $(GT) y) $(TD x - y)) |
| * $(TR $(TD x $(LT)= y) $(TD +0.0)) |
| * ) |
| */ |
| real fdim(real x, real y) @safe pure nothrow @nogc |
| { |
| return (x < y) ? +0.0 : x - y; |
| } |
| |
| /// |
| @safe pure nothrow @nogc unittest |
| { |
| import std.math.traits : isNaN; |
| |
| assert(fdim(2.0, 0.0) == 2.0); |
| assert(fdim(-2.0, 0.0) == 0.0); |
| assert(fdim(real.infinity, 2.0) == real.infinity); |
| assert(isNaN(fdim(real.nan, 2.0))); |
| assert(isNaN(fdim(2.0, real.nan))); |
| assert(isNaN(fdim(real.nan, real.nan))); |
| } |
| |
| /** |
| * Returns the larger of `x` and `y`. |
| * |
| * If one of the arguments is a `NaN`, the other is returned. |
| * |
| * See_Also: $(REF max, std,algorithm,comparison) is faster because it does not perform the `isNaN` test. |
| */ |
| F fmax(F)(const F x, const F y) @safe pure nothrow @nogc |
| if (__traits(isFloating, F)) |
| { |
| import std.math.traits : isNaN; |
| |
| // Do the more predictable test first. Generates 0 branches with ldc and 1 branch with gdc. |
| // See https://godbolt.org/z/erxrW9 |
| if (isNaN(x)) return y; |
| return y > x ? y : x; |
| } |
| |
| /// |
| @safe pure nothrow @nogc unittest |
| { |
| import std.meta : AliasSeq; |
| static foreach (F; AliasSeq!(float, double, real)) |
| { |
| assert(fmax(F(0.0), F(2.0)) == 2.0); |
| assert(fmax(F(-2.0), 0.0) == F(0.0)); |
| assert(fmax(F.infinity, F(2.0)) == F.infinity); |
| assert(fmax(F.nan, F(2.0)) == F(2.0)); |
| assert(fmax(F(2.0), F.nan) == F(2.0)); |
| } |
| } |
| |
| /** |
| * Returns the smaller of `x` and `y`. |
| * |
| * If one of the arguments is a `NaN`, the other is returned. |
| * |
| * See_Also: $(REF min, std,algorithm,comparison) is faster because it does not perform the `isNaN` test. |
| */ |
| F fmin(F)(const F x, const F y) @safe pure nothrow @nogc |
| if (__traits(isFloating, F)) |
| { |
| import std.math.traits : isNaN; |
| |
| // Do the more predictable test first. Generates 0 branches with ldc and 1 branch with gdc. |
| // See https://godbolt.org/z/erxrW9 |
| if (isNaN(x)) return y; |
| return y < x ? y : x; |
| } |
| |
| /// |
| @safe pure nothrow @nogc unittest |
| { |
| import std.meta : AliasSeq; |
| static foreach (F; AliasSeq!(float, double, real)) |
| { |
| assert(fmin(F(0.0), F(2.0)) == 0.0); |
| assert(fmin(F(-2.0), F(0.0)) == -2.0); |
| assert(fmin(F.infinity, F(2.0)) == 2.0); |
| assert(fmin(F.nan, F(2.0)) == 2.0); |
| assert(fmin(F(2.0), F.nan) == 2.0); |
| } |
| } |
| |
| /************************************** |
| * Returns (x * y) + z, rounding only once according to the |
| * current rounding mode. |
| * |
| * BUGS: Not currently implemented - rounds twice. |
| */ |
| pragma(inline, true) |
| real fma(real x, real y, real z) @safe pure nothrow @nogc { return (x * y) + z; } |
| |
| /// |
| @safe pure nothrow @nogc unittest |
| { |
| assert(fma(0.0, 2.0, 2.0) == 2.0); |
| assert(fma(2.0, 2.0, 2.0) == 6.0); |
| assert(fma(real.infinity, 2.0, 2.0) == real.infinity); |
| assert(fma(real.nan, 2.0, 2.0) is real.nan); |
| assert(fma(2.0, 2.0, real.nan) is real.nan); |
| } |
| |
| /************************************** |
| * To what precision is x equal to y? |
| * |
| * Returns: the number of mantissa bits which are equal in x and y. |
| * eg, 0x1.F8p+60 and 0x1.F1p+60 are equal to 5 bits of precision. |
| * |
| * $(TABLE_SV |
| * $(TR $(TH x) $(TH y) $(TH feqrel(x, y))) |
| * $(TR $(TD x) $(TD x) $(TD real.mant_dig)) |
| * $(TR $(TD x) $(TD $(GT)= 2*x) $(TD 0)) |
| * $(TR $(TD x) $(TD $(LT)= x/2) $(TD 0)) |
| * $(TR $(TD $(NAN)) $(TD any) $(TD 0)) |
| * $(TR $(TD any) $(TD $(NAN)) $(TD 0)) |
| * ) |
| */ |
| int feqrel(X)(const X x, const X y) @trusted pure nothrow @nogc |
| if (isFloatingPoint!(X)) |
| { |
| import std.math.traits : floatTraits, RealFormat; |
| import core.math : fabs; |
| |
| /* Public Domain. Author: Don Clugston, 18 Aug 2005. |
| */ |
| alias F = floatTraits!(X); |
| static if (F.realFormat == RealFormat.ieeeSingle |
| || F.realFormat == RealFormat.ieeeDouble |
| || F.realFormat == RealFormat.ieeeExtended |
| || F.realFormat == RealFormat.ieeeExtended53 |
| || F.realFormat == RealFormat.ieeeQuadruple) |
| { |
| if (x == y) |
| return X.mant_dig; // ensure diff != 0, cope with INF. |
| |
| Unqual!X diff = fabs(x - y); |
| |
| ushort *pa = cast(ushort *)(&x); |
| ushort *pb = cast(ushort *)(&y); |
| ushort *pd = cast(ushort *)(&diff); |
| |
| |
| // The difference in abs(exponent) between x or y and abs(x-y) |
| // is equal to the number of significand bits of x which are |
| // equal to y. If negative, x and y have different exponents. |
| // If positive, x and y are equal to 'bitsdiff' bits. |
| // AND with 0x7FFF to form the absolute value. |
| // To avoid out-by-1 errors, we subtract 1 so it rounds down |
| // if the exponents were different. This means 'bitsdiff' is |
| // always 1 lower than we want, except that if bitsdiff == 0, |
| // they could have 0 or 1 bits in common. |
| |
| int bitsdiff = ((( (pa[F.EXPPOS_SHORT] & F.EXPMASK) |
| + (pb[F.EXPPOS_SHORT] & F.EXPMASK) |
| - (1 << F.EXPSHIFT)) >> 1) |
| - (pd[F.EXPPOS_SHORT] & F.EXPMASK)) >> F.EXPSHIFT; |
| if ( (pd[F.EXPPOS_SHORT] & F.EXPMASK) == 0) |
| { // Difference is subnormal |
| // For subnormals, we need to add the number of zeros that |
| // lie at the start of diff's significand. |
| // We do this by multiplying by 2^^real.mant_dig |
| diff *= F.RECIP_EPSILON; |
| return bitsdiff + X.mant_dig - ((pd[F.EXPPOS_SHORT] & F.EXPMASK) >> F.EXPSHIFT); |
| } |
| |
| if (bitsdiff > 0) |
| return bitsdiff + 1; // add the 1 we subtracted before |
| |
| // Avoid out-by-1 errors when factor is almost 2. |
| if (bitsdiff == 0 |
| && ((pa[F.EXPPOS_SHORT] ^ pb[F.EXPPOS_SHORT]) & F.EXPMASK) == 0) |
| { |
| return 1; |
| } else return 0; |
| } |
| else |
| { |
| static assert(false, "Not implemented for this architecture"); |
| } |
| } |
| |
| /// |
| @safe pure unittest |
| { |
| assert(feqrel(2.0, 2.0) == 53); |
| assert(feqrel(2.0f, 2.0f) == 24); |
| assert(feqrel(2.0, double.nan) == 0); |
| |
| // Test that numbers are within n digits of each |
| // other by testing if feqrel > n * log2(10) |
| |
| // five digits |
| assert(feqrel(2.0, 2.00001) > 16); |
| // ten digits |
| assert(feqrel(2.0, 2.00000000001) > 33); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| void testFeqrel(F)() |
| { |
| // Exact equality |
| assert(feqrel(F.max, F.max) == F.mant_dig); |
| assert(feqrel!(F)(0.0, 0.0) == F.mant_dig); |
| assert(feqrel(F.infinity, F.infinity) == F.mant_dig); |
| |
| // a few bits away from exact equality |
| F w=1; |
| for (int i = 1; i < F.mant_dig - 1; ++i) |
| { |
| assert(feqrel!(F)(1.0 + w * F.epsilon, 1.0) == F.mant_dig-i); |
| assert(feqrel!(F)(1.0 - w * F.epsilon, 1.0) == F.mant_dig-i); |
| assert(feqrel!(F)(1.0, 1 + (w-1) * F.epsilon) == F.mant_dig - i + 1); |
| w*=2; |
| } |
| |
| assert(feqrel!(F)(1.5+F.epsilon, 1.5) == F.mant_dig-1); |
| assert(feqrel!(F)(1.5-F.epsilon, 1.5) == F.mant_dig-1); |
| assert(feqrel!(F)(1.5-F.epsilon, 1.5+F.epsilon) == F.mant_dig-2); |
| |
| |
| // Numbers that are close |
| assert(feqrel!(F)(0x1.Bp+84, 0x1.B8p+84) == 5); |
| assert(feqrel!(F)(0x1.8p+10, 0x1.Cp+10) == 2); |
| assert(feqrel!(F)(1.5 * (1 - F.epsilon), 1.0L) == 2); |
| assert(feqrel!(F)(1.5, 1.0) == 1); |
| assert(feqrel!(F)(2 * (1 - F.epsilon), 1.0L) == 1); |
| |
| // Factors of 2 |
| assert(feqrel(F.max, F.infinity) == 0); |
| assert(feqrel!(F)(2 * (1 - F.epsilon), 1.0L) == 1); |
| assert(feqrel!(F)(1.0, 2.0) == 0); |
| assert(feqrel!(F)(4.0, 1.0) == 0); |
| |
| // Extreme inequality |
| assert(feqrel(F.nan, F.nan) == 0); |
| assert(feqrel!(F)(0.0L, -F.nan) == 0); |
| assert(feqrel(F.nan, F.infinity) == 0); |
| assert(feqrel(F.infinity, -F.infinity) == 0); |
| assert(feqrel(F.max, -F.max) == 0); |
| |
| assert(feqrel(F.min_normal / 8, F.min_normal / 17) == 3); |
| |
| const F Const = 2; |
| immutable F Immutable = 2; |
| auto Compiles = feqrel(Const, Immutable); |
| } |
| |
| assert(feqrel(7.1824L, 7.1824L) == real.mant_dig); |
| |
| testFeqrel!(real)(); |
| testFeqrel!(double)(); |
| testFeqrel!(float)(); |
| } |
| |
| /** |
| Computes whether a values is approximately equal to a reference value, |
| admitting a maximum relative difference, and a maximum absolute difference. |
| |
| Warning: |
| This template is considered out-dated. It will be removed from |
| Phobos in 2.106.0. Please use $(LREF isClose) instead. To achieve |
| a similar behaviour to `approxEqual(a, b)` use |
| `isClose(a, b, 1e-2, 1e-5)`. In case of comparing to 0.0, |
| `isClose(a, b, 0.0, eps)` should be used, where `eps` |
| represents the accepted deviation from 0.0." |
| |
| Params: |
| value = Value to compare. |
| reference = Reference value. |
| maxRelDiff = Maximum allowable difference relative to `reference`. |
| Setting to 0.0 disables this check. Defaults to `1e-2`. |
| maxAbsDiff = Maximum absolute difference. This is mainly usefull |
| for comparing values to zero. Setting to 0.0 disables this check. |
| Defaults to `1e-5`. |
| |
| Returns: |
| `true` if `value` is approximately equal to `reference` under |
| either criterium. It is sufficient, when `value ` satisfies |
| one of the two criteria. |
| |
| If one item is a range, and the other is a single value, then |
| the result is the logical and-ing of calling `approxEqual` on |
| each element of the ranged item against the single item. If |
| both items are ranges, then `approxEqual` returns `true` if |
| and only if the ranges have the same number of elements and if |
| `approxEqual` evaluates to `true` for each pair of elements. |
| |
| See_Also: |
| Use $(LREF feqrel) to get the number of equal bits in the mantissa. |
| */ |
| deprecated("approxEqual will be removed in 2.106.0. Please use isClose instead.") |
| bool approxEqual(T, U, V)(T value, U reference, V maxRelDiff = 1e-2, V maxAbsDiff = 1e-5) |
| { |
| import core.math : fabs; |
| import std.range.primitives : empty, front, isInputRange, popFront; |
| static if (isInputRange!T) |
| { |
| static if (isInputRange!U) |
| { |
| // Two ranges |
| for (;; value.popFront(), reference.popFront()) |
| { |
| if (value.empty) return reference.empty; |
| if (reference.empty) return value.empty; |
| if (!approxEqual(value.front, reference.front, maxRelDiff, maxAbsDiff)) |
| return false; |
| } |
| } |
| else static if (isIntegral!U) |
| { |
| // convert reference to real |
| return approxEqual(value, real(reference), maxRelDiff, maxAbsDiff); |
| } |
| else |
| { |
| // value is range, reference is number |
| for (; !value.empty; value.popFront()) |
| { |
| if (!approxEqual(value.front, reference, maxRelDiff, maxAbsDiff)) |
| return false; |
| } |
| return true; |
| } |
| } |
| else |
| { |
| static if (isInputRange!U) |
| { |
| // value is number, reference is range |
| for (; !reference.empty; reference.popFront()) |
| { |
| if (!approxEqual(value, reference.front, maxRelDiff, maxAbsDiff)) |
| return false; |
| } |
| return true; |
| } |
| else static if (isIntegral!T || isIntegral!U) |
| { |
| // convert both value and reference to real |
| return approxEqual(real(value), real(reference), maxRelDiff, maxAbsDiff); |
| } |
| else |
| { |
| // two numbers |
| //static assert(is(T : real) && is(U : real)); |
| if (reference == 0) |
| { |
| return fabs(value) <= maxAbsDiff; |
| } |
| static if (is(typeof(value.infinity)) && is(typeof(reference.infinity))) |
| { |
| if (value == value.infinity && reference == reference.infinity || |
| value == -value.infinity && reference == -reference.infinity) return true; |
| } |
| return fabs((value - reference) / reference) <= maxRelDiff |
| || maxAbsDiff != 0 && fabs(value - reference) <= maxAbsDiff; |
| } |
| } |
| } |
| |
| deprecated @safe pure nothrow unittest |
| { |
| assert(approxEqual(1.0, 1.0099)); |
| assert(!approxEqual(1.0, 1.011)); |
| assert(approxEqual(0.00001, 0.0)); |
| assert(!approxEqual(0.00002, 0.0)); |
| |
| assert(approxEqual(3.0, [3, 3.01, 2.99])); // several reference values is strange |
| assert(approxEqual([3, 3.01, 2.99], 3.0)); // better |
| |
| float[] arr1 = [ 1.0, 2.0, 3.0 ]; |
| double[] arr2 = [ 1.001, 1.999, 3 ]; |
| assert(approxEqual(arr1, arr2)); |
| } |
| |
| deprecated @safe pure nothrow unittest |
| { |
| // relative comparison depends on reference, make sure proper |
| // side is used when comparing range to single value. Based on |
| // https://issues.dlang.org/show_bug.cgi?id=15763 |
| auto a = [2e-3 - 1e-5]; |
| auto b = 2e-3 + 1e-5; |
| assert(a[0].approxEqual(b)); |
| assert(!b.approxEqual(a[0])); |
| assert(a.approxEqual(b)); |
| assert(!b.approxEqual(a)); |
| } |
| |
| deprecated @safe pure nothrow @nogc unittest |
| { |
| assert(!approxEqual(0.0,1e-15,1e-9,0.0)); |
| assert(approxEqual(0.0,1e-15,1e-9,1e-9)); |
| assert(!approxEqual(1.0,3.0,0.0,1.0)); |
| |
| assert(approxEqual(1.00000000099,1.0,1e-9,0.0)); |
| assert(!approxEqual(1.0000000011,1.0,1e-9,0.0)); |
| } |
| |
| deprecated @safe pure nothrow @nogc unittest |
| { |
| // maybe unintuitive behavior |
| assert(approxEqual(1000.0,1010.0)); |
| assert(approxEqual(9_090_000_000.0,9_000_000_000.0)); |
| assert(approxEqual(0.0,1e30,1.0)); |
| assert(approxEqual(0.00001,1e-30)); |
| assert(!approxEqual(-1e-30,1e-30,1e-2,0.0)); |
| } |
| |
| deprecated @safe pure nothrow @nogc unittest |
| { |
| int a = 10; |
| assert(approxEqual(10, a)); |
| |
| assert(!approxEqual(3, 0)); |
| assert(approxEqual(3, 3)); |
| assert(approxEqual(3.0, 3)); |
| assert(approxEqual(3, 3.0)); |
| |
| assert(approxEqual(0.0,0.0)); |
| assert(approxEqual(-0.0,0.0)); |
| assert(approxEqual(0.0f,0.0)); |
| } |
| |
| deprecated @safe pure nothrow @nogc unittest |
| { |
| real num = real.infinity; |
| assert(num == real.infinity); |
| assert(approxEqual(num, real.infinity)); |
| num = -real.infinity; |
| assert(num == -real.infinity); |
| assert(approxEqual(num, -real.infinity)); |
| |
| assert(!approxEqual(1,real.nan)); |
| assert(!approxEqual(real.nan,real.max)); |
| assert(!approxEqual(real.nan,real.nan)); |
| } |
| |
| deprecated @safe pure nothrow unittest |
| { |
| assert(!approxEqual([1.0,2.0,3.0],[1.0,2.0])); |
| assert(!approxEqual([1.0,2.0],[1.0,2.0,3.0])); |
| |
| assert(approxEqual!(real[],real[])([],[])); |
| assert(approxEqual(cast(real[])[],cast(real[])[])); |
| } |
| |
| |
| /** |
| Computes whether two values are approximately equal, admitting a maximum |
| relative difference, and a maximum absolute difference. |
| |
| Params: |
| lhs = First item to compare. |
| rhs = Second item to compare. |
| maxRelDiff = Maximum allowable relative difference. |
| Setting to 0.0 disables this check. Default depends on the type of |
| `lhs` and `rhs`: It is approximately half the number of decimal digits of |
| precision of the smaller type. |
| maxAbsDiff = Maximum absolute difference. This is mainly usefull |
| for comparing values to zero. Setting to 0.0 disables this check. |
| Defaults to `0.0`. |
| |
| Returns: |
| `true` if the two items are approximately equal under either criterium. |
| It is sufficient, when `value ` satisfies one of the two criteria. |
| |
| If one item is a range, and the other is a single value, then |
| the result is the logical and-ing of calling `isClose` on |
| each element of the ranged item against the single item. If |
| both items are ranges, then `isClose` returns `true` if |
| and only if the ranges have the same number of elements and if |
| `isClose` evaluates to `true` for each pair of elements. |
| |
| See_Also: |
| Use $(LREF feqrel) to get the number of equal bits in the mantissa. |
| */ |
| bool isClose(T, U, V = CommonType!(FloatingPointBaseType!T,FloatingPointBaseType!U)) |
| (T lhs, U rhs, V maxRelDiff = CommonDefaultFor!(T,U), V maxAbsDiff = 0.0) |
| { |
| import std.range.primitives : empty, front, isInputRange, popFront; |
| import std.complex : Complex; |
| static if (isInputRange!T) |
| { |
| static if (isInputRange!U) |
| { |
| // Two ranges |
| for (;; lhs.popFront(), rhs.popFront()) |
| { |
| if (lhs.empty) return rhs.empty; |
| if (rhs.empty) return lhs.empty; |
| if (!isClose(lhs.front, rhs.front, maxRelDiff, maxAbsDiff)) |
| return false; |
| } |
| } |
| else |
| { |
| // lhs is range, rhs is number |
| for (; !lhs.empty; lhs.popFront()) |
| { |
| if (!isClose(lhs.front, rhs, maxRelDiff, maxAbsDiff)) |
| return false; |
| } |
| return true; |
| } |
| } |
| else static if (isInputRange!U) |
| { |
| // lhs is number, rhs is range |
| for (; !rhs.empty; rhs.popFront()) |
| { |
| if (!isClose(lhs, rhs.front, maxRelDiff, maxAbsDiff)) |
| return false; |
| } |
| return true; |
| } |
| else static if (is(T TE == Complex!TE)) |
| { |
| static if (is(U UE == Complex!UE)) |
| { |
| // Two complex numbers |
| return isClose(lhs.re, rhs.re, maxRelDiff, maxAbsDiff) |
| && isClose(lhs.im, rhs.im, maxRelDiff, maxAbsDiff); |
| } |
| else |
| { |
| // lhs is complex, rhs is number |
| return isClose(lhs.re, rhs, maxRelDiff, maxAbsDiff) |
| && isClose(lhs.im, 0.0, maxRelDiff, maxAbsDiff); |
| } |
| } |
| else static if (is(U UE == Complex!UE)) |
| { |
| // lhs is number, rhs is complex |
| return isClose(lhs, rhs.re, maxRelDiff, maxAbsDiff) |
| && isClose(0.0, rhs.im, maxRelDiff, maxAbsDiff); |
| } |
| else |
| { |
| // two numbers |
| if (lhs == rhs) return true; |
| |
| static if (is(typeof(lhs.infinity))) |
| if (lhs == lhs.infinity || lhs == -lhs.infinity) |
| return false; |
| static if (is(typeof(rhs.infinity))) |
| if (rhs == rhs.infinity || rhs == -rhs.infinity) |
| return false; |
| |
| import std.math.algebraic : abs; |
| |
| auto diff = abs(lhs - rhs); |
| |
| return diff <= maxRelDiff*abs(lhs) |
| || diff <= maxRelDiff*abs(rhs) |
| || diff <= maxAbsDiff; |
| } |
| } |
| |
| /// |
| @safe pure nothrow @nogc unittest |
| { |
| assert(isClose(1.0,0.999_999_999)); |
| assert(isClose(0.001, 0.000_999_999_999)); |
| assert(isClose(1_000_000_000.0,999_999_999.0)); |
| |
| assert(isClose(17.123_456_789, 17.123_456_78)); |
| assert(!isClose(17.123_456_789, 17.123_45)); |
| |
| // use explicit 3rd parameter for less (or more) accuracy |
| assert(isClose(17.123_456_789, 17.123_45, 1e-6)); |
| assert(!isClose(17.123_456_789, 17.123_45, 1e-7)); |
| |
| // use 4th parameter when comparing close to zero |
| assert(!isClose(1e-100, 0.0)); |
| assert(isClose(1e-100, 0.0, 0.0, 1e-90)); |
| assert(!isClose(1e-10, -1e-10)); |
| assert(isClose(1e-10, -1e-10, 0.0, 1e-9)); |
| assert(!isClose(1e-300, 1e-298)); |
| assert(isClose(1e-300, 1e-298, 0.0, 1e-200)); |
| |
| // different default limits for different floating point types |
| assert(isClose(1.0f, 0.999_99f)); |
| assert(!isClose(1.0, 0.999_99)); |
| static if (real.sizeof > double.sizeof) |
| assert(!isClose(1.0L, 0.999_999_999L)); |
| } |
| |
| /// |
| @safe pure nothrow unittest |
| { |
| assert(isClose([1.0, 2.0, 3.0], [0.999_999_999, 2.000_000_001, 3.0])); |
| assert(!isClose([1.0, 2.0], [0.999_999_999, 2.000_000_001, 3.0])); |
| assert(!isClose([1.0, 2.0, 3.0], [0.999_999_999, 2.000_000_001])); |
| |
| assert(isClose([2.0, 1.999_999_999, 2.000_000_001], 2.0)); |
| assert(isClose(2.0, [2.0, 1.999_999_999, 2.000_000_001])); |
| } |
| |
| @safe pure nothrow unittest |
| { |
| assert(!isClose([1.0, 2.0, 3.0], [0.999_999_999, 3.0, 3.0])); |
| assert(!isClose([2.0, 1.999_999, 2.000_000_001], 2.0)); |
| assert(!isClose(2.0, [2.0, 1.999_999_999, 2.000_000_999])); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| immutable a = 1.00001f; |
| const b = 1.000019; |
| assert(isClose(a,b)); |
| |
| assert(isClose(1.00001f,1.000019f)); |
| assert(isClose(1.00001f,1.000019)); |
| assert(isClose(1.00001,1.000019f)); |
| assert(!isClose(1.00001,1.000019)); |
| |
| real a1 = 1e-300L; |
| real a2 = a1.nextUp; |
| assert(isClose(a1,a2)); |
| } |
| |
| @safe pure nothrow unittest |
| { |
| float[] arr1 = [ 1.0, 2.0, 3.0 ]; |
| double[] arr2 = [ 1.00001, 1.99999, 3 ]; |
| assert(isClose(arr1, arr2)); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| assert(!isClose(1000.0,1010.0)); |
| assert(!isClose(9_090_000_000.0,9_000_000_000.0)); |
| assert(isClose(0.0,1e30,1.0)); |
| assert(!isClose(0.00001,1e-30)); |
| assert(!isClose(-1e-30,1e-30,1e-2,0.0)); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| assert(!isClose(3, 0)); |
| assert(isClose(3, 3)); |
| assert(isClose(3.0, 3)); |
| assert(isClose(3, 3.0)); |
| |
| assert(isClose(0.0,0.0)); |
| assert(isClose(-0.0,0.0)); |
| assert(isClose(0.0f,0.0)); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| real num = real.infinity; |
| assert(num == real.infinity); |
| assert(isClose(num, real.infinity)); |
| num = -real.infinity; |
| assert(num == -real.infinity); |
| assert(isClose(num, -real.infinity)); |
| |
| assert(!isClose(1,real.nan)); |
| assert(!isClose(real.nan,real.max)); |
| assert(!isClose(real.nan,real.nan)); |
| |
| assert(!isClose(-double.infinity, 1)); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| assert(isClose!(real[],real[],real)([],[])); |
| assert(isClose(cast(real[])[],cast(real[])[])); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| import std.conv : to; |
| |
| float f = 31.79f; |
| double d = 31.79; |
| double f2d = f.to!double; |
| |
| assert(isClose(f,f2d)); |
| assert(!isClose(d,f2d)); |
| } |
| |
| @safe pure nothrow @nogc unittest |
| { |
| import std.conv : to; |
| |
| double d = 31.79; |
| float f = d.to!float; |
| double f2d = f.to!double; |
| |
| assert(isClose(f,f2d)); |
| assert(!isClose(d,f2d)); |
| assert(isClose(d,f2d,1e-4)); |
| } |
| |
| package(std.math) template CommonDefaultFor(T,U) |
| { |
| import std.algorithm.comparison : min; |
| |
| alias baseT = FloatingPointBaseType!T; |
| alias baseU = FloatingPointBaseType!U; |
| |
| enum CommonType!(baseT, baseU) CommonDefaultFor = 10.0L ^^ -((min(baseT.dig, baseU.dig) + 1) / 2 + 1); |
| } |
| |
| private template FloatingPointBaseType(T) |
| { |
| import std.range.primitives : ElementType; |
| static if (isFloatingPoint!T) |
| { |
| alias FloatingPointBaseType = Unqual!T; |
| } |
| else static if (isFloatingPoint!(ElementType!(Unqual!T))) |
| { |
| alias FloatingPointBaseType = Unqual!(ElementType!(Unqual!T)); |
| } |
| else |
| { |
| alias FloatingPointBaseType = real; |
| } |
| } |
| |
| /*********************************** |
| * Defines a total order on all floating-point numbers. |
| * |
| * The order is defined as follows: |
| * $(UL |
| * $(LI All numbers in [-$(INFIN), +$(INFIN)] are ordered |
| * the same way as by built-in comparison, with the exception of |
| * -0.0, which is less than +0.0;) |
| * $(LI If the sign bit is set (that is, it's 'negative'), $(NAN) is less |
| * than any number; if the sign bit is not set (it is 'positive'), |
| * $(NAN) is greater than any number;) |
| * $(LI $(NAN)s of the same sign are ordered by the payload ('negative' |
| * ones - in reverse order).) |
| * ) |
| * |
| * Returns: |
| * negative value if `x` precedes `y` in the order specified above; |
| * 0 if `x` and `y` are identical, and positive value otherwise. |
| * |
| * See_Also: |
| * $(MYREF isIdentical) |
| * Standards: Conforms to IEEE 754-2008 |
| */ |
| int cmp(T)(const(T) x, const(T) y) @nogc @trusted pure nothrow |
| if (isFloatingPoint!T) |
| { |
| import std.math.traits : floatTraits, RealFormat; |
| |
| alias F = floatTraits!T; |
| |
| static if (F.realFormat == RealFormat.ieeeSingle |
| || F.realFormat == RealFormat.ieeeDouble) |
| { |
| static if (T.sizeof == 4) |
| alias UInt = uint; |
| else |
| alias UInt = ulong; |
| |
| union Repainter |
| { |
| T number; |
| UInt bits; |
| } |
| |
| enum msb = ~(UInt.max >>> 1); |
| |
| import std.typecons : Tuple; |
| Tuple!(Repainter, Repainter) vars = void; |
| vars[0].number = x; |
| vars[1].number = y; |
| |
| foreach (ref var; vars) |
| if (var.bits & msb) |
| var.bits = ~var.bits; |
| else |
| var.bits |= msb; |
| |
| if (vars[0].bits < vars[1].bits) |
| return -1; |
| else if (vars[0].bits > vars[1].bits) |
| return 1; |
| else |
| return 0; |
| } |
| else static if (F.realFormat == RealFormat.ieeeExtended53 |
| || F.realFormat == RealFormat.ieeeExtended |
| || F.realFormat == RealFormat.ieeeQuadruple) |
| { |
| static if (F.realFormat == RealFormat.ieeeQuadruple) |
| alias RemT = ulong; |
| else |
| alias RemT = ushort; |
| |
| struct Bits |
| { |
| ulong bulk; |
| RemT rem; |
| } |
| |
| union Repainter |
| { |
| T number; |
| Bits bits; |
| ubyte[T.sizeof] bytes; |
| } |
| |
| import std.typecons : Tuple; |
| Tuple!(Repainter, Repainter) vars = void; |
| vars[0].number = x; |
| vars[1].number = y; |
| |
| foreach (ref var; vars) |
| if (var.bytes[F.SIGNPOS_BYTE] & 0x80) |
| { |
| var.bits.bulk = ~var.bits.bulk; |
| var.bits.rem = cast(typeof(var.bits.rem))(-1 - var.bits.rem); // ~var.bits.rem |
| } |
| else |
| { |
| var.bytes[F.SIGNPOS_BYTE] |= 0x80; |
| } |
| |
| version (LittleEndian) |
| { |
| if (vars[0].bits.rem < vars[1].bits.rem) |
| return -1; |
| else if (vars[0].bits.rem > vars[1].bits.rem) |
| return 1; |
| else if (vars[0].bits.bulk < vars[1].bits.bulk) |
| return -1; |
| else if (vars[0].bits.bulk > vars[1].bits.bulk) |
| return 1; |
| else |
| return 0; |
| } |
| else |
| { |
| if (vars[0].bits.bulk < vars[1].bits.bulk) |
| return -1; |
| else if (vars[0].bits.bulk > vars[1].bits.bulk) |
| return 1; |
| else if (vars[0].bits.rem < vars[1].bits.rem) |
| return -1; |
| else if (vars[0].bits.rem > vars[1].bits.rem) |
| return 1; |
| else |
| return 0; |
| } |
| } |
| else |
| { |
| // IBM Extended doubledouble does not follow the general |
| // sign-exponent-significand layout, so has to be handled generically |
| |
| import std.math.traits : signbit, isNaN; |
| |
| const int xSign = signbit(x), |
| ySign = signbit(y); |
| |
| if (xSign == 1 && ySign == 1) |
| return cmp(-y, -x); |
| else if (xSign == 1) |
| return -1; |
| else if (ySign == 1) |
| return 1; |
| else if (x < y) |
| return -1; |
| else if (x == y) |
| return 0; |
| else if (x > y) |
| return 1; |
| else if (isNaN(x) && !isNaN(y)) |
| return 1; |
| else if (isNaN(y) && !isNaN(x)) |
| return -1; |
| else if (getNaNPayload(x) < getNaNPayload(y)) |
| return -1; |
| else if (getNaNPayload(x) > getNaNPayload(y)) |
| return 1; |
| else |
| return 0; |
| } |
| } |
| |
| /// Most numbers are ordered naturally. |
| @safe unittest |
| { |
| assert(cmp(-double.infinity, -double.max) < 0); |
| assert(cmp(-double.max, -100.0) < 0); |
| assert(cmp(-100.0, -0.5) < 0); |
| assert(cmp(-0.5, 0.0) < 0); |
| assert(cmp(0.0, 0.5) < 0); |
| assert(cmp(0.5, 100.0) < 0); |
| assert(cmp(100.0, double.max) < 0); |
| assert(cmp(double.max, double.infinity) < 0); |
| |
| assert(cmp(1.0, 1.0) == 0); |
| } |
| |
| /// Positive and negative zeroes are distinct. |
| @safe unittest |
| { |
| assert(cmp(-0.0, +0.0) < 0); |
| assert(cmp(+0.0, -0.0) > 0); |
| } |
| |
| /// Depending on the sign, $(NAN)s go to either end of the spectrum. |
| @safe unittest |
| { |
| assert(cmp(-double.nan, -double.infinity) < 0); |
| assert(cmp(double.infinity, double.nan) < 0); |
| assert(cmp(-double.nan, double.nan) < 0); |
| } |
| |
| /// $(NAN)s of the same sign are ordered by the payload. |
| @safe unittest |
| { |
| assert(cmp(NaN(10), NaN(20)) < 0); |
| assert(cmp(-NaN(20), -NaN(10)) < 0); |
| } |
| |
| @safe unittest |
| { |
| import std.meta : AliasSeq; |
| static foreach (T; AliasSeq!(float, double, real)) |
| {{ |
| T[] values = [-cast(T) NaN(20), -cast(T) NaN(10), -T.nan, -T.infinity, |
| -T.max, -T.max / 2, T(-16.0), T(-1.0).nextDown, |
| T(-1.0), T(-1.0).nextUp, |
| T(-0.5), -T.min_normal, (-T.min_normal).nextUp, |
| -2 * T.min_normal * T.epsilon, |
| -T.min_normal * T.epsilon, |
| T(-0.0), T(0.0), |
| T.min_normal * T.epsilon, |
| 2 * T.min_normal * T.epsilon, |
| T.min_normal.nextDown, T.min_normal, T(0.5), |
| T(1.0).nextDown, T(1.0), |
| T(1.0).nextUp, T(16.0), T.max / 2, T.max, |
| T.infinity, T.nan, cast(T) NaN(10), cast(T) NaN(20)]; |
| |
| foreach (i, x; values) |
| { |
| foreach (y; values[i + 1 .. $]) |
| { |
| assert(cmp(x, y) < 0); |
| assert(cmp(y, x) > 0); |
| } |
| assert(cmp(x, x) == 0); |
| } |
| }} |
| } |
| |
| package(std): // not yet public |
| |
| struct FloatingPointBitpattern(T) |
| if (isFloatingPoint!T) |
| { |
| static if (T.mant_dig <= 64) |
| { |
| ulong mantissa; |
| } |
| else |
| { |
| ulong mantissa_lsb; |
| ulong mantissa_msb; |
| } |
| |
| int exponent; |
| bool negative; |
| } |
| |
| FloatingPointBitpattern!T extractBitpattern(T)(const(T) value) @trusted |
| if (isFloatingPoint!T) |
| { |
| import std.math.traits : floatTraits, RealFormat; |
| |
| T val = value; |
| FloatingPointBitpattern!T ret; |
| |
| alias F = floatTraits!T; |
| static if (F.realFormat == RealFormat.ieeeExtended) |
| { |
| if (__ctfe) |
| { |
| import core.math : fabs, ldexp; |
| import std.math.rounding : floor; |
| import std.math.traits : isInfinity, isNaN, signbit; |
| import std.math.exponential : log2; |
| |
| if (isNaN(val) || isInfinity(val)) |
| ret.exponent = 32767; |
| else if (fabs(val) < real.min_normal) |
| ret.exponent = 0; |
| else if (fabs(val) >= nextUp(real.max / 2)) |
| ret.exponent = 32766; |
| else |
| ret.exponent = cast(int) (val.fabs.log2.floor() + 16383); |
| |
| if (ret.exponent == 32767) |
| { |
| // NaN or infinity |
| ret.mantissa = isNaN(val) ? ((1L << 63) - 1) : 0; |
| } |
| else |
| { |
| auto delta = 16382 + 64 // bias + bits of ulong |
| - (ret.exponent == 0 ? 1 : ret.exponent); // -1 in case of subnormals |
| val = ldexp(val, delta); // val *= 2^^delta |
| |
| ulong tmp = cast(ulong) fabs(val); |
| if (ret.exponent != 32767 && ret.exponent > 0 && tmp <= ulong.max / 2) |
| { |
| // correction, due to log2(val) being rounded up: |
| ret.exponent--; |
| val *= 2; |
| tmp = cast(ulong) fabs(val); |
| } |
| |
| ret.mantissa = tmp & long.max; |
| } |
| |
| ret.negative = (signbit(val) == 1); |
| } |
| else |
| { |
| ushort* vs = cast(ushort*) &val; |
| ret.mantissa = (cast(ulong*) vs)[0] & long.max; |
| ret.exponent = vs[4] & short.max; |
| ret.negative = (vs[4] >> 15) & 1; |
| } |
| } |
| else |
| { |
| static if (F.realFormat == RealFormat.ieeeSingle) |
| { |
| ulong ival = *cast(uint*) &val; |
| } |
| else static if (F.realFormat == RealFormat.ieeeDouble) |
| { |
| ulong ival = *cast(ulong*) &val; |
| } |
| else |
| { |
| static assert(false, "Floating point type `" ~ F.realFormat ~ "` not supported."); |
| } |
| |
| import std.math.exponential : log2; |
| enum log2_max_exp = cast(int) log2(T(T.max_exp)); |
| |
| ret.mantissa = ival & ((1L << (T.mant_dig - 1)) - 1); |
| ret.exponent = (ival >> (T.mant_dig - 1)) & ((1L << (log2_max_exp + 1)) - 1); |
| ret.negative = (ival >> (T.mant_dig + log2_max_exp)) & 1; |
| } |
| |
| // add leading 1 for normalized values and correct exponent for denormalied values |
| if (ret.exponent != 0 && ret.exponent != 2 * T.max_exp - 1) |
| ret.mantissa |= 1L << (T.mant_dig - 1); |
| else if (ret.exponent == 0) |
| ret.exponent = 1; |
| |
| ret.exponent -= T.max_exp - 1; |
| |
| return ret; |
| } |
| |
| @safe pure unittest |
| { |
| float f = 1.0f; |
| auto bp = extractBitpattern(f); |
| assert(bp.mantissa == 0x80_0000); |
| assert(bp.exponent == 0); |
| assert(bp.negative == false); |
| |
| f = float.max; |
| bp = extractBitpattern(f); |
| assert(bp.mantissa == 0xff_ffff); |
| assert(bp.exponent == 127); |
| assert(bp.negative == false); |
| |
| f = -1.5432e-17f; |
| bp = extractBitpattern(f); |
| assert(bp.mantissa == 0x8e_55c8); |
| assert(bp.exponent == -56); |
| assert(bp.negative == true); |
| |
| // using double literal due to https://issues.dlang.org/show_bug.cgi?id=20361 |
| f = 2.3822073893521890206e-44; |
| bp = extractBitpattern(f); |
| assert(bp.mantissa == 0x00_0011); |
| assert(bp.exponent == -126); |
| assert(bp.negative == false); |
| |
| f = -float.infinity; |
| bp = extractBitpattern(f); |
| assert(bp.mantissa == 0); |
| assert(bp.exponent == 128); |
| assert(bp.negative == true); |
| |
| f = float.nan; |
| bp = extractBitpattern(f); |
| assert(bp.mantissa != 0); // we don't guarantee payloads |
| assert(bp.exponent == 128); |
| assert(bp.negative == false); |
| } |
| |
| @safe pure unittest |
| { |
| double d = 1.0; |
| auto bp = extractBitpattern(d); |
| assert(bp.mantissa == 0x10_0000_0000_0000L); |
| assert(bp.exponent == 0); |
| assert(bp.negative == false); |
| |
| d = double.max; |
| bp = extractBitpattern(d); |
| assert(bp.mantissa == 0x1f_ffff_ffff_ffffL); |
| assert(bp.exponent == 1023); |
| assert(bp.negative == false); |
| |
| d = -1.5432e-222; |
| bp = extractBitpattern(d); |
| assert(bp.mantissa == 0x11_d9b6_a401_3b04L); |
| assert(bp.exponent == -737); |
| assert(bp.negative == true); |
| |
| d = 0.0.nextUp; |
| bp = extractBitpattern(d); |
| assert(bp.mantissa == 0x00_0000_0000_0001L); |
| assert(bp.exponent == -1022); |
| assert(bp.negative == false); |
| |
| d = -double.infinity; |
| bp = extractBitpattern(d); |
| assert(bp.mantissa == 0); |
| assert(bp.exponent == 1024); |
| assert(bp.negative == true); |
| |
| d = double.nan; |
| bp = extractBitpattern(d); |
| assert(bp.mantissa != 0); // we don't guarantee payloads |
| assert(bp.exponent == 1024); |
| assert(bp.negative == false); |
| } |
| |
| @safe pure unittest |
| { |
| import std.math.traits : floatTraits, RealFormat; |
| |
| alias F = floatTraits!real; |
| static if (F.realFormat == RealFormat.ieeeExtended) |
| { |
| real r = 1.0L; |
| auto bp = extractBitpattern(r); |
| assert(bp.mantissa == 0x8000_0000_0000_0000L); |
| assert(bp.exponent == 0); |
| assert(bp.negative == false); |
| |
| r = real.max; |
| bp = extractBitpattern(r); |
| assert(bp.mantissa == 0xffff_ffff_ffff_ffffL); |
| assert(bp.exponent == 16383); |
| assert(bp.negative == false); |
| |
| r = -1.5432e-3333L; |
| bp = extractBitpattern(r); |
| assert(bp.mantissa == 0xc768_a2c7_a616_cc22L); |
| assert(bp.exponent == -11072); |
| assert(bp.negative == true); |
| |
| r = 0.0L.nextUp; |
| bp = extractBitpattern(r); |
| assert(bp.mantissa == 0x0000_0000_0000_0001L); |
| assert(bp.exponent == -16382); |
| assert(bp.negative == false); |
| |
| r = -float.infinity; |
| bp = extractBitpattern(r); |
| assert(bp.mantissa == 0); |
| assert(bp.exponent == 16384); |
| assert(bp.negative == true); |
| |
| r = float.nan; |
| bp = extractBitpattern(r); |
| assert(bp.mantissa != 0); // we don't guarantee payloads |
| assert(bp.exponent == 16384); |
| assert(bp.negative == false); |
| |
| r = nextDown(0x1p+16383L); |
| bp = extractBitpattern(r); |
| assert(bp.mantissa == 0xffff_ffff_ffff_ffffL); |
| assert(bp.exponent == 16382); |
| assert(bp.negative == false); |
| } |
| } |
| |
| @safe pure unittest |
| { |
| import std.math.traits : floatTraits, RealFormat; |
| import std.math.exponential : log2; |
| |
| alias F = floatTraits!real; |
| |
| static if (F.realFormat == RealFormat.ieeeExtended) |
| { |
| // log2 is broken for x87-reals on some computers in CTFE |
| // the following test excludes these computers from the test |
| // (https://issues.dlang.org/show_bug.cgi?id=21757) |
| enum test = cast(int) log2(3.05e2312L); |
| static if (test == 7681) |
| { |
| enum r1 = 1.0L; |
| enum bp1 = extractBitpattern(r1); |
| static assert(bp1.mantissa == 0x8000_0000_0000_0000L); |
| static assert(bp1.exponent == 0); |
| static assert(bp1.negative == false); |
| |
| enum r2 = real.max; |
| enum bp2 = extractBitpattern(r2); |
| static assert(bp2.mantissa == 0xffff_ffff_ffff_ffffL); |
| static assert(bp2.exponent == 16383); |
| static assert(bp2.negative == false); |
| |
| enum r3 = -1.5432e-3333L; |
| enum bp3 = extractBitpattern(r3); |
| static assert(bp3.mantissa == 0xc768_a2c7_a616_cc22L); |
| static assert(bp3.exponent == -11072); |
| static assert(bp3.negative == true); |
| |
| enum r4 = 0.0L.nextUp; |
| enum bp4 = extractBitpattern(r4); |
| static assert(bp4.mantissa == 0x0000_0000_0000_0001L); |
| static assert(bp4.exponent == -16382); |
| static assert(bp4.negative == false); |
| |
| enum r5 = -real.infinity; |
| enum bp5 = extractBitpattern(r5); |
| static assert(bp5.mantissa == 0); |
| static assert(bp5.exponent == 16384); |
| static assert(bp5.negative == true); |
| |
| enum r6 = real.nan; |
| enum bp6 = extractBitpattern(r6); |
| static assert(bp6.mantissa != 0); // we don't guarantee payloads |
| static assert(bp6.exponent == 16384); |
| static assert(bp6.negative == false); |
| |
| enum r7 = nextDown(0x1p+16383L); |
| enum bp7 = extractBitpattern(r7); |
| static assert(bp7.mantissa == 0xffff_ffff_ffff_ffffL); |
| static assert(bp7.exponent == 16382); |
| static assert(bp7.negative == false); |
| } |
| } |
| } |
