libphobos/src/std/math/traits.d

   1 // Written in the D programming language.
   2
   3 /**
   4 This is a submodule of $(MREF std, math).
   5
   6 It contains several functions for introspection on numerical values.
   7
   8 Copyright: Copyright The D Language Foundation 2000 - 2011.
   9 License:   $(HTTP www.boost.org/LICENSE_1_0.txt, Boost License 1.0).
  10 Authors:   $(HTTP digitalmars.com, Walter Bright), Don Clugston,
  11            Conversion of CEPHES math library to D by Iain Buclaw and David Nadlinger
  12 Source: $(PHOBOSSRC std/math/traits.d)
  13
  14 Macros:
  15     NAN = $(RED NAN)
  16     PLUSMN = &plusmn;
  17     INFIN = &infin;
  18  */
  19
  20 module std.math.traits;
  21
  22 import std.traits : isFloatingPoint, isIntegral, isNumeric, isSigned;
  23
  24
  25 /*********************************
  26  * Determines if $(D_PARAM x) is NaN.
  27  * Params:
  28  *  x = a floating point number.
  29  * Returns:
  30  *  `true` if $(D_PARAM x) is Nan.
  31  */
  32 bool isNaN(X)(X x) @nogc @trusted pure nothrow
  33 if (isFloatingPoint!(X))
  34 {
  35     version (all)
  36     {
  37         return x != x;
  38     }
  39     else
  40     {
  41         /*
  42         Code kept for historical context. At least on Intel, the simple test
  43         x != x uses one dedicated instruction (ucomiss/ucomisd) that runs in one
  44         cycle. Code for 80- and 128-bits is larger but still smaller than the
  45         integrals-based solutions below. Future revisions may enable the code
  46         below conditionally depending on hardware.
  47         */
  48         alias F = floatTraits!(X);
  49         static if (F.realFormat == RealFormat.ieeeSingle)
  50         {
  51             const uint p = *cast(uint *)&x;
  52             // Sign bit (MSB) is irrelevant so mask it out.
  53             // Next 8 bits should be all set.
  54             // At least one bit among the least significant 23 bits should be set.
  55             return (p & 0x7FFF_FFFF) > 0x7F80_0000;
  56         }
  57         else static if (F.realFormat == RealFormat.ieeeDouble)
  58         {
  59             const ulong  p = *cast(ulong *)&x;
  60             // Sign bit (MSB) is irrelevant so mask it out.
  61             // Next 11 bits should be all set.
  62             // At least one bit among the least significant 52 bits should be set.
  63             return (p & 0x7FFF_FFFF_FFFF_FFFF) > 0x7FF0_0000_0000_0000;
  64         }
  65         else static if (F.realFormat == RealFormat.ieeeExtended ||
  66                         F.realFormat == RealFormat.ieeeExtended53)
  67         {
  68             const ushort e = F.EXPMASK & (cast(ushort *)&x)[F.EXPPOS_SHORT];
  69             const ulong ps = *cast(ulong *)&x;
  70             return e == F.EXPMASK &&
  71                 ps & 0x7FFF_FFFF_FFFF_FFFF; // not infinity
  72         }
  73         else static if (F.realFormat == RealFormat.ieeeQuadruple)
  74         {
  75             const ushort e = F.EXPMASK & (cast(ushort *)&x)[F.EXPPOS_SHORT];
  76             const ulong psLsb = (cast(ulong *)&x)[MANTISSA_LSB];
  77             const ulong psMsb = (cast(ulong *)&x)[MANTISSA_MSB];
  78             return e == F.EXPMASK &&
  79                 (psLsb | (psMsb& 0x0000_FFFF_FFFF_FFFF)) != 0;
  80         }
  81         else
  82         {
  83             return x != x;
  84         }
  85     }
  86 }
  87
  88 ///
  89 @safe pure nothrow @nogc unittest
  90 {
  91     assert( isNaN(float.init));
  92     assert( isNaN(-double.init));
  93     assert( isNaN(real.nan));
  94     assert( isNaN(-real.nan));
  95     assert(!isNaN(cast(float) 53.6));
  96     assert(!isNaN(cast(real)-53.6));
  97 }
  98
  99 @safe pure nothrow @nogc unittest
 100 {
 101     import std.meta : AliasSeq;
 102
 103     static foreach (T; AliasSeq!(float, double, real))
 104     {{
 105         // CTFE-able tests
 106         assert(isNaN(T.init));
 107         assert(isNaN(-T.init));
 108         assert(isNaN(T.nan));
 109         assert(isNaN(-T.nan));
 110         assert(!isNaN(T.infinity));
 111         assert(!isNaN(-T.infinity));
 112         assert(!isNaN(cast(T) 53.6));
 113         assert(!isNaN(cast(T)-53.6));
 114
 115         // Runtime tests
 116         shared T f;
 117         f = T.init;
 118         assert(isNaN(f));
 119         assert(isNaN(-f));
 120         f = T.nan;
 121         assert(isNaN(f));
 122         assert(isNaN(-f));
 123         f = T.infinity;
 124         assert(!isNaN(f));
 125         assert(!isNaN(-f));
 126         f = cast(T) 53.6;
 127         assert(!isNaN(f));
 128         assert(!isNaN(-f));
 129     }}
 130 }
 131
 132 /*********************************
 133  * Determines if $(D_PARAM x) is finite.
 134  * Params:
 135  *  x = a floating point number.
 136  * Returns:
 137  *  `true` if $(D_PARAM x) is finite.
 138  */
 139 bool isFinite(X)(X x) @trusted pure nothrow @nogc
 140 {
 141     import std.math.traits : floatTraits, RealFormat;
 142
 143     static if (__traits(isFloating, X))
 144         if (__ctfe)
 145             return x == x && x != X.infinity && x != -X.infinity;
 146     alias F = floatTraits!(X);
 147     ushort* pe = cast(ushort *)&x;
 148     return (pe[F.EXPPOS_SHORT] & F.EXPMASK) != F.EXPMASK;
 149 }
 150
 151 ///
 152 @safe pure nothrow @nogc unittest
 153 {
 154     assert( isFinite(1.23f));
 155     assert( isFinite(float.max));
 156     assert( isFinite(float.min_normal));
 157     assert(!isFinite(float.nan));
 158     assert(!isFinite(float.infinity));
 159 }
 160
 161 @safe pure nothrow @nogc unittest
 162 {
 163     assert(isFinite(1.23));
 164     assert(isFinite(double.max));
 165     assert(isFinite(double.min_normal));
 166     assert(!isFinite(double.nan));
 167     assert(!isFinite(double.infinity));
 168
 169     assert(isFinite(1.23L));
 170     assert(isFinite(real.max));
 171     assert(isFinite(real.min_normal));
 172     assert(!isFinite(real.nan));
 173     assert(!isFinite(real.infinity));
 174
 175     //CTFE
 176     static assert(isFinite(1.23));
 177     static assert(isFinite(double.max));
 178     static assert(isFinite(double.min_normal));
 179     static assert(!isFinite(double.nan));
 180     static assert(!isFinite(double.infinity));
 181
 182     static assert(isFinite(1.23L));
 183     static assert(isFinite(real.max));
 184     static assert(isFinite(real.min_normal));
 185     static assert(!isFinite(real.nan));
 186     static assert(!isFinite(real.infinity));
 187 }
 188
 189
 190 /*********************************
 191  * Determines if $(D_PARAM x) is normalized.
 192  *
 193  * A normalized number must not be zero, subnormal, infinite nor $(NAN).
 194  *
 195  * Params:
 196  *  x = a floating point number.
 197  * Returns:
 198  *  `true` if $(D_PARAM x) is normalized.
 199  */
 200
 201 /* Need one for each format because subnormal floats might
 202  * be converted to normal reals.
 203  */
 204 bool isNormal(X)(X x) @trusted pure nothrow @nogc
 205 {
 206     import std.math.traits : floatTraits, RealFormat;
 207
 208     static if (__traits(isFloating, X))
 209         if (__ctfe)
 210             return (x <= -X.min_normal && x != -X.infinity) || (x >= X.min_normal && x != X.infinity);
 211     alias F = floatTraits!(X);
 212     ushort e = F.EXPMASK & (cast(ushort *)&x)[F.EXPPOS_SHORT];
 213     return (e != F.EXPMASK && e != 0);
 214 }
 215
 216 ///
 217 @safe pure nothrow @nogc unittest
 218 {
 219     float f = 3;
 220     double d = 500;
 221     real e = 10e+48;
 222
 223     assert(isNormal(f));
 224     assert(isNormal(d));
 225     assert(isNormal(e));
 226     f = d = e = 0;
 227     assert(!isNormal(f));
 228     assert(!isNormal(d));
 229     assert(!isNormal(e));
 230     assert(!isNormal(real.infinity));
 231     assert(isNormal(-real.max));
 232     assert(!isNormal(real.min_normal/4));
 233
 234 }
 235
 236 @safe pure nothrow @nogc unittest
 237 {
 238     // CTFE
 239     enum float f = 3;
 240     enum double d = 500;
 241     enum real e = 10e+48;
 242
 243     static assert(isNormal(f));
 244     static assert(isNormal(d));
 245     static assert(isNormal(e));
 246
 247     static assert(!isNormal(0.0f));
 248     static assert(!isNormal(0.0));
 249     static assert(!isNormal(0.0L));
 250     static assert(!isNormal(real.infinity));
 251     static assert(isNormal(-real.max));
 252     static assert(!isNormal(real.min_normal/4));
 253 }
 254
 255 /*********************************
 256  * Determines if $(D_PARAM x) is subnormal.
 257  *
 258  * Subnormals (also known as "denormal number"), have a 0 exponent
 259  * and a 0 most significant mantissa bit.
 260  *
 261  * Params:
 262  *  x = a floating point number.
 263  * Returns:
 264  *  `true` if $(D_PARAM x) is a denormal number.
 265  */
 266 bool isSubnormal(X)(X x) @trusted pure nothrow @nogc
 267 {
 268     import std.math.traits : floatTraits, RealFormat, MANTISSA_MSB, MANTISSA_LSB;
 269
 270     static if (__traits(isFloating, X))
 271         if (__ctfe)
 272             return -X.min_normal < x && x < X.min_normal;
 273     /*
 274         Need one for each format because subnormal floats might
 275         be converted to normal reals.
 276     */
 277     alias F = floatTraits!(X);
 278     static if (F.realFormat == RealFormat.ieeeSingle)
 279     {
 280         uint *p = cast(uint *)&x;
 281         return (*p & F.EXPMASK_INT) == 0 && *p & F.MANTISSAMASK_INT;
 282     }
 283     else static if (F.realFormat == RealFormat.ieeeDouble)
 284     {
 285         uint *p = cast(uint *)&x;
 286         return (p[MANTISSA_MSB] & F.EXPMASK_INT) == 0
 287             && (p[MANTISSA_LSB] || p[MANTISSA_MSB] & F.MANTISSAMASK_INT);
 288     }
 289     else static if (F.realFormat == RealFormat.ieeeQuadruple)
 290     {
 291         ushort e = F.EXPMASK & (cast(ushort *)&x)[F.EXPPOS_SHORT];
 292         long*   ps = cast(long *)&x;
 293         return (e == 0 &&
 294           ((ps[MANTISSA_LSB]|(ps[MANTISSA_MSB]& 0x0000_FFFF_FFFF_FFFF)) != 0));
 295     }
 296     else static if (F.realFormat == RealFormat.ieeeExtended ||
 297                     F.realFormat == RealFormat.ieeeExtended53)
 298     {
 299         ushort* pe = cast(ushort *)&x;
 300         long*   ps = cast(long *)&x;
 301
 302         return (pe[F.EXPPOS_SHORT] & F.EXPMASK) == 0 && *ps > 0;
 303     }
 304     else
 305     {
 306         static assert(false, "Not implemented for this architecture");
 307     }
 308 }
 309
 310 ///
 311 @safe pure nothrow @nogc unittest
 312 {
 313     import std.meta : AliasSeq;
 314
 315     static foreach (T; AliasSeq!(float, double, real))
 316     {{
 317         T f;
 318         for (f = 1.0; !isSubnormal(f); f /= 2)
 319             assert(f != 0);
 320     }}
 321 }
 322
 323 @safe pure nothrow @nogc unittest
 324 {
 325     static bool subnormalTest(T)()
 326     {
 327         T f;
 328         for (f = 1.0; !isSubnormal(f); f /= 2)
 329             if (f == 0)
 330                 return false;
 331         return true;
 332     }
 333     static assert(subnormalTest!float());
 334     static assert(subnormalTest!double());
 335     static assert(subnormalTest!real());
 336 }
 337
 338 /*********************************
 339  * Determines if $(D_PARAM x) is $(PLUSMN)$(INFIN).
 340  * Params:
 341  *  x = a floating point number.
 342  * Returns:
 343  *  `true` if $(D_PARAM x) is $(PLUSMN)$(INFIN).
 344  */
 345 bool isInfinity(X)(X x) @nogc @trusted pure nothrow
 346 if (isFloatingPoint!(X))
 347 {
 348     import std.math.traits : floatTraits, RealFormat, MANTISSA_MSB, MANTISSA_LSB;
 349
 350     alias F = floatTraits!(X);
 351     static if (F.realFormat == RealFormat.ieeeSingle)
 352     {
 353         return ((*cast(uint *)&x) & 0x7FFF_FFFF) == 0x7F80_0000;
 354     }
 355     else static if (F.realFormat == RealFormat.ieeeDouble)
 356     {
 357         return ((*cast(ulong *)&x) & 0x7FFF_FFFF_FFFF_FFFF)
 358             == 0x7FF0_0000_0000_0000;
 359     }
 360     else static if (F.realFormat == RealFormat.ieeeExtended ||
 361                     F.realFormat == RealFormat.ieeeExtended53)
 362     {
 363         const ushort e = cast(ushort)(F.EXPMASK & (cast(ushort *)&x)[F.EXPPOS_SHORT]);
 364         const ulong ps = *cast(ulong *)&x;
 365
 366         // On Motorola 68K, infinity can have hidden bit = 1 or 0. On x86, it is always 1.
 367         return e == F.EXPMASK && (ps & 0x7FFF_FFFF_FFFF_FFFF) == 0;
 368     }
 369     else static if (F.realFormat == RealFormat.ieeeQuadruple)
 370     {
 371         const long psLsb = (cast(long *)&x)[MANTISSA_LSB];
 372         const long psMsb = (cast(long *)&x)[MANTISSA_MSB];
 373         return (psLsb == 0)
 374             && (psMsb & 0x7FFF_FFFF_FFFF_FFFF) == 0x7FFF_0000_0000_0000;
 375     }
 376     else
 377     {
 378         return (x < -X.max) || (X.max < x);
 379     }
 380 }
 381
 382 ///
 383 @nogc @safe pure nothrow unittest
 384 {
 385     assert(!isInfinity(float.init));
 386     assert(!isInfinity(-float.init));
 387     assert(!isInfinity(float.nan));
 388     assert(!isInfinity(-float.nan));
 389     assert(isInfinity(float.infinity));
 390     assert(isInfinity(-float.infinity));
 391     assert(isInfinity(-1.0f / 0.0f));
 392 }
 393
 394 @safe pure nothrow @nogc unittest
 395 {
 396     // CTFE-able tests
 397     assert(!isInfinity(double.init));
 398     assert(!isInfinity(-double.init));
 399     assert(!isInfinity(double.nan));
 400     assert(!isInfinity(-double.nan));
 401     assert(isInfinity(double.infinity));
 402     assert(isInfinity(-double.infinity));
 403     assert(isInfinity(-1.0 / 0.0));
 404
 405     assert(!isInfinity(real.init));
 406     assert(!isInfinity(-real.init));
 407     assert(!isInfinity(real.nan));
 408     assert(!isInfinity(-real.nan));
 409     assert(isInfinity(real.infinity));
 410     assert(isInfinity(-real.infinity));
 411     assert(isInfinity(-1.0L / 0.0L));
 412
 413     // Runtime tests
 414     shared float f;
 415     f = float.init;
 416     assert(!isInfinity(f));
 417     assert(!isInfinity(-f));
 418     f = float.nan;
 419     assert(!isInfinity(f));
 420     assert(!isInfinity(-f));
 421     f = float.infinity;
 422     assert(isInfinity(f));
 423     assert(isInfinity(-f));
 424     f = (-1.0f / 0.0f);
 425     assert(isInfinity(f));
 426
 427     shared double d;
 428     d = double.init;
 429     assert(!isInfinity(d));
 430     assert(!isInfinity(-d));
 431     d = double.nan;
 432     assert(!isInfinity(d));
 433     assert(!isInfinity(-d));
 434     d = double.infinity;
 435     assert(isInfinity(d));
 436     assert(isInfinity(-d));
 437     d = (-1.0 / 0.0);
 438     assert(isInfinity(d));
 439
 440     shared real e;
 441     e = real.init;
 442     assert(!isInfinity(e));
 443     assert(!isInfinity(-e));
 444     e = real.nan;
 445     assert(!isInfinity(e));
 446     assert(!isInfinity(-e));
 447     e = real.infinity;
 448     assert(isInfinity(e));
 449     assert(isInfinity(-e));
 450     e = (-1.0L / 0.0L);
 451     assert(isInfinity(e));
 452 }
 453
 454 @nogc @safe pure nothrow unittest
 455 {
 456     import std.meta : AliasSeq;
 457     static bool foo(T)(inout T x) { return isInfinity(x); }
 458     foreach (T; AliasSeq!(float, double, real))
 459     {
 460         assert(!foo(T(3.14f)));
 461         assert(foo(T.infinity));
 462     }
 463 }
 464
 465 /*********************************
 466  * Is the binary representation of x identical to y?
 467  */
 468 bool isIdentical(real x, real y) @trusted pure nothrow @nogc
 469 {
 470     if (__ctfe)
 471     {
 472         if (x !is y) return false;
 473         if (x == x) return true; // If not NaN `is` implies identical representation.
 474         static if (double.mant_dig != real.mant_dig)
 475         {
 476             // Works because we are in CTFE and there is no way in CTFE to set more
 477             // bits of NaN payload than can fit in a double, and since 2.087
 478             // changed real.init to be non-signaling I *think* there is no way in
 479             // CTFE for a real to be a signaling NaN unless real and double have
 480             // the same representation so real's bits can be manipulated directly.
 481             double d1 = x, d2 = y;
 482         }
 483         else
 484         {
 485             // Alias to avoid converting signaling to quiet.
 486             alias d1 = x;
 487             alias d2 = y;
 488         }
 489         return *cast(long*) &d1 == *cast(long*) &d2;
 490     }
 491     // We're doing a bitwise comparison so the endianness is irrelevant.
 492     long*   pxs = cast(long *)&x;
 493     long*   pys = cast(long *)&y;
 494     alias F = floatTraits!(real);
 495     static if (F.realFormat == RealFormat.ieeeDouble)
 496     {
 497         return pxs[0] == pys[0];
 498     }
 499     else static if (F.realFormat == RealFormat.ieeeQuadruple)
 500     {
 501         return pxs[0] == pys[0] && pxs[1] == pys[1];
 502     }
 503     else static if (F.realFormat == RealFormat.ieeeExtended)
 504     {
 505         ushort* pxe = cast(ushort *)&x;
 506         ushort* pye = cast(ushort *)&y;
 507         return pxe[4] == pye[4] && pxs[0] == pys[0];
 508     }
 509     else
 510     {
 511         assert(0, "isIdentical not implemented");
 512     }
 513 }
 514 ///
 515 @safe @nogc pure nothrow unittest
 516 {
 517     // We're forcing the CTFE to run by assigning the result of the function to an enum
 518     enum test1 = isIdentical(1.0,1.0);
 519     enum test2 = isIdentical(real.nan,real.nan);
 520     enum test3 = isIdentical(real.infinity, real.infinity);
 521     enum test4 = isIdentical(real.infinity, real.infinity);
 522     enum test5 = isIdentical(0.0, 0.0);
 523
 524     assert(test1);
 525     assert(test2);
 526     assert(test3);
 527     assert(test4);
 528     assert(test5);
 529
 530     enum test6 = !isIdentical(0.0, -0.0);
 531     enum test7 = !isIdentical(real.nan, -real.nan);
 532     enum test8 = !isIdentical(real.infinity, -real.infinity);
 533
 534     assert(test6);
 535     assert(test7);
 536     assert(test8);
 537 }
 538
 539 @safe @nogc pure nothrow unittest
 540 {
 541     assert( !isIdentical(1.2,1.3));
 542     assert( isIdentical(0.0, 0.0));
 543     assert( isIdentical(1.0, 1.0));
 544     assert( isIdentical(real.infinity, real.infinity));
 545     assert( isIdentical(-real.infinity, -real.infinity));
 546     assert( isIdentical(real.nan, real.nan));
 547
 548     assert(!isIdentical(0.0, -0.0));
 549     assert(!isIdentical(real.nan, -real.nan));
 550     assert(!isIdentical(real.infinity, -real.infinity));
 551 }
 552 /*********************************
 553  * Return 1 if sign bit of e is set, 0 if not.
 554  */
 555 int signbit(X)(X x) @nogc @trusted pure nothrow
 556 {
 557     import std.math.traits : floatTraits, RealFormat;
 558
 559     if (__ctfe)
 560     {
 561         double dval = cast(double) x; // Precision can increase or decrease but sign won't change (even NaN).
 562         return 0 > *cast(long*) &dval;
 563     }
 564
 565     alias F = floatTraits!(X);
 566     return ((cast(ubyte *)&x)[F.SIGNPOS_BYTE] & 0x80) != 0;
 567 }
 568
 569 ///
 570 @nogc @safe pure nothrow unittest
 571 {
 572     assert(!signbit(float.nan));
 573     assert(signbit(-float.nan));
 574     assert(!signbit(168.1234f));
 575     assert(signbit(-168.1234f));
 576     assert(!signbit(0.0f));
 577     assert(signbit(-0.0f));
 578     assert(signbit(-float.max));
 579     assert(!signbit(float.max));
 580
 581     assert(!signbit(double.nan));
 582     assert(signbit(-double.nan));
 583     assert(!signbit(168.1234));
 584     assert(signbit(-168.1234));
 585     assert(!signbit(0.0));
 586     assert(signbit(-0.0));
 587     assert(signbit(-double.max));
 588     assert(!signbit(double.max));
 589
 590     assert(!signbit(real.nan));
 591     assert(signbit(-real.nan));
 592     assert(!signbit(168.1234L));
 593     assert(signbit(-168.1234L));
 594     assert(!signbit(0.0L));
 595     assert(signbit(-0.0L));
 596     assert(signbit(-real.max));
 597     assert(!signbit(real.max));
 598 }
 599
 600 @nogc @safe pure nothrow unittest
 601 {
 602     // CTFE
 603     static assert(!signbit(float.nan));
 604     static assert(signbit(-float.nan));
 605     static assert(!signbit(168.1234f));
 606     static assert(signbit(-168.1234f));
 607     static assert(!signbit(0.0f));
 608     static assert(signbit(-0.0f));
 609     static assert(signbit(-float.max));
 610     static assert(!signbit(float.max));
 611
 612     static assert(!signbit(double.nan));
 613     static assert(signbit(-double.nan));
 614     static assert(!signbit(168.1234));
 615     static assert(signbit(-168.1234));
 616     static assert(!signbit(0.0));
 617     static assert(signbit(-0.0));
 618     static assert(signbit(-double.max));
 619     static assert(!signbit(double.max));
 620
 621     static assert(!signbit(real.nan));
 622     static assert(signbit(-real.nan));
 623     static assert(!signbit(168.1234L));
 624     static assert(signbit(-168.1234L));
 625     static assert(!signbit(0.0L));
 626     static assert(signbit(-0.0L));
 627     static assert(signbit(-real.max));
 628     static assert(!signbit(real.max));
 629 }
 630
 631 /**
 632 Params:
 633     to = the numeric value to use
 634     from = the sign value to use
 635 Returns:
 636     a value composed of to with from's sign bit.
 637  */
 638 R copysign(R, X)(R to, X from) @trusted pure nothrow @nogc
 639 if (isFloatingPoint!(R) && isFloatingPoint!(X))
 640 {
 641     import std.math.traits : floatTraits, RealFormat;
 642
 643     if (__ctfe)
 644     {
 645         return signbit(to) == signbit(from) ? to : -to;
 646     }
 647     ubyte* pto   = cast(ubyte *)&to;
 648     const ubyte* pfrom = cast(ubyte *)&from;
 649
 650     alias T = floatTraits!(R);
 651     alias F = floatTraits!(X);
 652     pto[T.SIGNPOS_BYTE] &= 0x7F;
 653     pto[T.SIGNPOS_BYTE] |= pfrom[F.SIGNPOS_BYTE] & 0x80;
 654     return to;
 655 }
 656
 657 /// ditto
 658 R copysign(R, X)(X to, R from) @trusted pure nothrow @nogc
 659 if (isIntegral!(X) && isFloatingPoint!(R))
 660 {
 661     return copysign(cast(R) to, from);
 662 }
 663
 664 ///
 665 @safe pure nothrow @nogc unittest
 666 {
 667     assert(copysign(1.0, 1.0) == 1.0);
 668     assert(copysign(1.0, -0.0) == -1.0);
 669     assert(copysign(1UL, -1.0) == -1.0);
 670     assert(copysign(-1.0, -1.0) == -1.0);
 671
 672     assert(copysign(real.infinity, -1.0) == -real.infinity);
 673     assert(copysign(real.nan, 1.0) is real.nan);
 674     assert(copysign(-real.nan, 1.0) is real.nan);
 675     assert(copysign(real.nan, -1.0) is -real.nan);
 676 }
 677
 678 @safe pure nothrow @nogc unittest
 679 {
 680     import std.meta : AliasSeq;
 681
 682     static foreach (X; AliasSeq!(float, double, real, int, long))
 683     {
 684         static foreach (Y; AliasSeq!(float, double, real))
 685         {{
 686             X x = 21;
 687             Y y = 23.8;
 688             Y e = void;
 689
 690             e = copysign(x, y);
 691             assert(e == 21.0);
 692
 693             e = copysign(-x, y);
 694             assert(e == 21.0);
 695
 696             e = copysign(x, -y);
 697             assert(e == -21.0);
 698
 699             e = copysign(-x, -y);
 700             assert(e == -21.0);
 701
 702             static if (isFloatingPoint!X)
 703             {
 704                 e = copysign(X.nan, y);
 705                 assert(isNaN(e) && !signbit(e));
 706
 707                 e = copysign(X.nan, -y);
 708                 assert(isNaN(e) && signbit(e));
 709             }
 710         }}
 711     }
 712     // CTFE
 713     static foreach (X; AliasSeq!(float, double, real, int, long))
 714     {
 715         static foreach (Y; AliasSeq!(float, double, real))
 716         {{
 717             enum X x = 21;
 718             enum Y y = 23.8;
 719
 720             assert(21.0 == copysign(x, y));
 721             assert(21.0 == copysign(-x, y));
 722             assert(-21.0 == copysign(x, -y));
 723             assert(-21.0 == copysign(-x, -y));
 724
 725             static if (isFloatingPoint!X)
 726             {
 727                 static assert(isNaN(copysign(X.nan, y)) && !signbit(copysign(X.nan, y)));
 728                 assert(isNaN(copysign(X.nan, -y)) && signbit(copysign(X.nan, -y)));
 729             }
 730         }}
 731     }
 732 }
 733
 734 /*********************************
 735 Returns `-1` if $(D x < 0), `x` if $(D x == 0), `1` if
 736 $(D x > 0), and $(NAN) if x==$(NAN).
 737  */
 738 F sgn(F)(F x) @safe pure nothrow @nogc
 739 if (isFloatingPoint!F || isIntegral!F)
 740 {
 741     // @@@TODO@@@: make this faster
 742     return x > 0 ? 1 : x < 0 ? -1 : x;
 743 }
 744
 745 ///
 746 @safe pure nothrow @nogc unittest
 747 {
 748     assert(sgn(168.1234) == 1);
 749     assert(sgn(-168.1234) == -1);
 750     assert(sgn(0.0) == 0);
 751     assert(sgn(-0.0) == 0);
 752 }
 753
 754 /**
 755 Check whether a number is an integer power of two.
 756
 757 Note that only positive numbers can be integer powers of two. This
 758 function always return `false` if `x` is negative or zero.
 759
 760 Params:
 761     x = the number to test
 762
 763 Returns:
 764     `true` if `x` is an integer power of two.
 765 */
 766 bool isPowerOf2(X)(const X x) pure @safe nothrow @nogc
 767 if (isNumeric!X)
 768 {
 769     import std.math.exponential : frexp;
 770
 771     static if (isFloatingPoint!X)
 772     {
 773         int exp;
 774         const X sig = frexp(x, exp);
 775
 776         return (exp != int.min) && (sig is cast(X) 0.5L);
 777     }
 778     else
 779     {
 780         static if (isSigned!X)
 781         {
 782             auto y = cast(typeof(x + 0))x;
 783             return y > 0 && !(y & (y - 1));
 784         }
 785         else
 786         {
 787             auto y = cast(typeof(x + 0u))x;
 788             return (y & -y) > (y - 1);
 789         }
 790     }
 791 }
 792 ///
 793 @safe unittest
 794 {
 795     import std.math.exponential : pow;
 796
 797     assert( isPowerOf2(1.0L));
 798     assert( isPowerOf2(2.0L));
 799     assert( isPowerOf2(0.5L));
 800     assert( isPowerOf2(pow(2.0L, 96)));
 801     assert( isPowerOf2(pow(2.0L, -77)));
 802
 803     assert(!isPowerOf2(-2.0L));
 804     assert(!isPowerOf2(-0.5L));
 805     assert(!isPowerOf2(0.0L));
 806     assert(!isPowerOf2(4.315));
 807     assert(!isPowerOf2(1.0L / 3.0L));
 808
 809     assert(!isPowerOf2(real.nan));
 810     assert(!isPowerOf2(real.infinity));
 811 }
 812 ///
 813 @safe unittest
 814 {
 815     assert( isPowerOf2(1));
 816     assert( isPowerOf2(2));
 817     assert( isPowerOf2(1uL << 63));
 818
 819     assert(!isPowerOf2(-4));
 820     assert(!isPowerOf2(0));
 821     assert(!isPowerOf2(1337u));
 822 }
 823
 824 @safe unittest
 825 {
 826     import std.math.exponential : pow;
 827     import std.meta : AliasSeq;
 828
 829     enum smallP2 = pow(2.0L, -62);
 830     enum bigP2 = pow(2.0L, 50);
 831     enum smallP7 = pow(7.0L, -35);
 832     enum bigP7 = pow(7.0L, 30);
 833
 834     static foreach (X; AliasSeq!(float, double, real))
 835     {{
 836         immutable min_sub = X.min_normal * X.epsilon;
 837
 838         foreach (x; [smallP2, min_sub, X.min_normal, .25L, 0.5L, 1.0L,
 839                               2.0L, 8.0L, pow(2.0L, X.max_exp - 1), bigP2])
 840         {
 841             assert( isPowerOf2(cast(X) x));
 842             assert(!isPowerOf2(cast(X)-x));
 843         }
 844
 845         foreach (x; [0.0L, 3 * min_sub, smallP7, 0.1L, 1337.0L, bigP7, X.max, real.nan, real.infinity])
 846         {
 847             assert(!isPowerOf2(cast(X) x));
 848             assert(!isPowerOf2(cast(X)-x));
 849         }
 850     }}
 851
 852     static foreach (X; AliasSeq!(byte, ubyte, short, ushort, int, uint, long, ulong))
 853     {{
 854         foreach (x; [1, 2, 4, 8, (X.max >>> 1) + 1])
 855         {
 856             assert( isPowerOf2(cast(X) x));
 857             static if (isSigned!X)
 858                 assert(!isPowerOf2(cast(X)-x));
 859         }
 860
 861         foreach (x; [0, 3, 5, 13, 77, X.min, X.max])
 862             assert(!isPowerOf2(cast(X) x));
 863     }}
 864
 865     // CTFE
 866     static foreach (X; AliasSeq!(float, double, real))
 867     {{
 868         enum min_sub = X.min_normal * X.epsilon;
 869
 870         static foreach (x; [smallP2, min_sub, X.min_normal, .25L, 0.5L, 1.0L,
 871                               2.0L, 8.0L, pow(2.0L, X.max_exp - 1), bigP2])
 872         {
 873             static assert( isPowerOf2(cast(X) x));
 874             static assert(!isPowerOf2(cast(X)-x));
 875         }
 876
 877         static foreach (x; [0.0L, 3 * min_sub, smallP7, 0.1L, 1337.0L, bigP7, X.max, real.nan, real.infinity])
 878         {
 879             static assert(!isPowerOf2(cast(X) x));
 880             static assert(!isPowerOf2(cast(X)-x));
 881         }
 882     }}
 883
 884     static foreach (X; AliasSeq!(byte, ubyte, short, ushort, int, uint, long, ulong))
 885     {{
 886         static foreach (x; [1, 2, 4, 8, (X.max >>> 1) + 1])
 887         {
 888             static assert( isPowerOf2(cast(X) x));
 889             static if (isSigned!X)
 890                 static assert(!isPowerOf2(cast(X)-x));
 891         }
 892
 893         static foreach (x; [0, 3, 5, 13, 77, X.min, X.max])
 894             static assert(!isPowerOf2(cast(X) x));
 895     }}
 896 }
 897
 898 // Underlying format exposed through floatTraits
 899 enum RealFormat
 900 {
 901     ieeeHalf,
 902     ieeeSingle,
 903     ieeeDouble,
 904     ieeeExtended,   // x87 80-bit real
 905     ieeeExtended53, // x87 real rounded to precision of double.
 906     ibmExtended,    // IBM 128-bit extended
 907     ieeeQuadruple,
 908 }
 909
 910 // Constants used for extracting the components of the representation.
 911 // They supplement the built-in floating point properties.
 912 template floatTraits(T)
 913 {
 914     import std.traits : Unqual;
 915
 916     // EXPMASK is a ushort mask to select the exponent portion (without sign)
 917     // EXPSHIFT is the number of bits the exponent is left-shifted by in its ushort
 918     // EXPBIAS is the exponent bias - 1 (exp == EXPBIAS yields ×2^-1).
 919     // EXPPOS_SHORT is the index of the exponent when represented as a ushort array.
 920     // SIGNPOS_BYTE is the index of the sign when represented as a ubyte array.
 921     // RECIP_EPSILON is the value such that (smallest_subnormal) * RECIP_EPSILON == T.min_normal
 922     enum Unqual!T RECIP_EPSILON = (1/T.epsilon);
 923     static if (T.mant_dig == 24)
 924     {
 925         // Single precision float
 926         enum ushort EXPMASK = 0x7F80;
 927         enum ushort EXPSHIFT = 7;
 928         enum ushort EXPBIAS = 0x3F00;
 929         enum uint EXPMASK_INT = 0x7F80_0000;
 930         enum uint MANTISSAMASK_INT = 0x007F_FFFF;
 931         enum realFormat = RealFormat.ieeeSingle;
 932         version (LittleEndian)
 933         {
 934             enum EXPPOS_SHORT = 1;
 935             enum SIGNPOS_BYTE = 3;
 936         }
 937         else
 938         {
 939             enum EXPPOS_SHORT = 0;
 940             enum SIGNPOS_BYTE = 0;
 941         }
 942     }
 943     else static if (T.mant_dig == 53)
 944     {
 945         static if (T.sizeof == 8)
 946         {
 947             // Double precision float, or real == double
 948             enum ushort EXPMASK = 0x7FF0;
 949             enum ushort EXPSHIFT = 4;
 950             enum ushort EXPBIAS = 0x3FE0;
 951             enum uint EXPMASK_INT = 0x7FF0_0000;
 952             enum uint MANTISSAMASK_INT = 0x000F_FFFF; // for the MSB only
 953             enum ulong MANTISSAMASK_LONG = 0x000F_FFFF_FFFF_FFFF;
 954             enum realFormat = RealFormat.ieeeDouble;
 955             version (LittleEndian)
 956             {
 957                 enum EXPPOS_SHORT = 3;
 958                 enum SIGNPOS_BYTE = 7;
 959             }
 960             else
 961             {
 962                 enum EXPPOS_SHORT = 0;
 963                 enum SIGNPOS_BYTE = 0;
 964             }
 965         }
 966         else static if (T.sizeof == 12)
 967         {
 968             // Intel extended real80 rounded to double
 969             enum ushort EXPMASK = 0x7FFF;
 970             enum ushort EXPSHIFT = 0;
 971             enum ushort EXPBIAS = 0x3FFE;
 972             enum realFormat = RealFormat.ieeeExtended53;
 973             version (LittleEndian)
 974             {
 975                 enum EXPPOS_SHORT = 4;
 976                 enum SIGNPOS_BYTE = 9;
 977             }
 978             else
 979             {
 980                 enum EXPPOS_SHORT = 0;
 981                 enum SIGNPOS_BYTE = 0;
 982             }
 983         }
 984         else
 985             static assert(false, "No traits support for " ~ T.stringof);
 986     }
 987     else static if (T.mant_dig == 64)
 988     {
 989         // Intel extended real80
 990         enum ushort EXPMASK = 0x7FFF;
 991         enum ushort EXPSHIFT = 0;
 992         enum ushort EXPBIAS = 0x3FFE;
 993         enum realFormat = RealFormat.ieeeExtended;
 994         version (LittleEndian)
 995         {
 996             enum EXPPOS_SHORT = 4;
 997             enum SIGNPOS_BYTE = 9;
 998         }
 999         else
1000         {
1001             enum EXPPOS_SHORT = 0;
1002             enum SIGNPOS_BYTE = 0;
1003         }
1004     }
1005     else static if (T.mant_dig == 113)
1006     {
1007         // Quadruple precision float
1008         enum ushort EXPMASK = 0x7FFF;
1009         enum ushort EXPSHIFT = 0;
1010         enum ushort EXPBIAS = 0x3FFE;
1011         enum realFormat = RealFormat.ieeeQuadruple;
1012         version (LittleEndian)
1013         {
1014             enum EXPPOS_SHORT = 7;
1015             enum SIGNPOS_BYTE = 15;
1016         }
1017         else
1018         {
1019             enum EXPPOS_SHORT = 0;
1020             enum SIGNPOS_BYTE = 0;
1021         }
1022     }
1023     else static if (T.mant_dig == 106)
1024     {
1025         // IBM Extended doubledouble
1026         enum ushort EXPMASK = 0x7FF0;
1027         enum ushort EXPSHIFT = 4;
1028         enum realFormat = RealFormat.ibmExtended;
1029
1030         // For IBM doubledouble the larger magnitude double comes first.
1031         // It's really a double[2] and arrays don't index differently
1032         // between little and big-endian targets.
1033         enum DOUBLEPAIR_MSB = 0;
1034         enum DOUBLEPAIR_LSB = 1;
1035
1036         // The exponent/sign byte is for most significant part.
1037         version (LittleEndian)
1038         {
1039             enum EXPPOS_SHORT = 3;
1040             enum SIGNPOS_BYTE = 7;
1041         }
1042         else
1043         {
1044             enum EXPPOS_SHORT = 0;
1045             enum SIGNPOS_BYTE = 0;
1046         }
1047     }
1048     else
1049         static assert(false, "No traits support for " ~ T.stringof);
1050 }
1051
1052 // These apply to all floating-point types
1053 version (LittleEndian)
1054 {
1055     enum MANTISSA_LSB = 0;
1056     enum MANTISSA_MSB = 1;
1057 }
1058 else
1059 {
1060     enum MANTISSA_LSB = 1;
1061     enum MANTISSA_MSB = 0;
1062 }