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 rounding floating point numbers. 7 8 Copyright: Copyright The D Language Foundation 2000 - 2011. 9 D implementations of floor, ceil, and lrint functions are based on the 10 CEPHES math library, which is Copyright (C) 2001 Stephen L. Moshier 11 $(LT)steve@moshier.net$(GT) and are incorporated herein by permission 12 of the author. The author reserves the right to distribute this 13 material elsewhere under different copying permissions. 14 These modifications are distributed here under the following terms: 15 License: $(HTTP www.boost.org/LICENSE_1_0.txt, Boost License 1.0). 16 Authors: $(HTTP digitalmars.com, Walter Bright), Don Clugston, 17 Conversion of CEPHES math library to D by Iain Buclaw and David Nadlinger 18 Source: $(PHOBOSSRC std/math/rounding.d) 19 */ 20 21 module std.math.rounding; 22 23 static import core.math; 24 static import core.stdc.math; 25 26 import std.traits : isFloatingPoint, isIntegral, Unqual; 27 28 version (D_InlineAsm_X86) version = InlineAsm_X86_Any; 29 version (D_InlineAsm_X86_64) version = InlineAsm_X86_Any; 30 31 version (InlineAsm_X86_Any) version = InlineAsm_X87; 32 version (InlineAsm_X87) 33 { 34 static assert(real.mant_dig == 64); 35 version (CRuntime_Microsoft) version = InlineAsm_X87_MSVC; 36 } 37 38 /************************************** 39 * Returns the value of x rounded upward to the next integer 40 * (toward positive infinity). 41 */ 42 real ceil(real x) @trusted pure nothrow @nogc 43 { 44 version (InlineAsm_X87_MSVC) 45 { 46 version (X86_64) 47 { 48 asm pure nothrow @nogc 49 { 50 naked ; 51 fld real ptr [RCX] ; 52 fstcw 8[RSP] ; 53 mov AL,9[RSP] ; 54 mov DL,AL ; 55 and AL,0xC3 ; 56 or AL,0x08 ; // round to +infinity 57 mov 9[RSP],AL ; 58 fldcw 8[RSP] ; 59 frndint ; 60 mov 9[RSP],DL ; 61 fldcw 8[RSP] ; 62 ret ; 63 } 64 } 65 else 66 { 67 short cw; 68 asm pure nothrow @nogc 69 { 70 fld x ; 71 fstcw cw ; 72 mov AL,byte ptr cw+1 ; 73 mov DL,AL ; 74 and AL,0xC3 ; 75 or AL,0x08 ; // round to +infinity 76 mov byte ptr cw+1,AL ; 77 fldcw cw ; 78 frndint ; 79 mov byte ptr cw+1,DL ; 80 fldcw cw ; 81 } 82 } 83 } 84 else 85 { 86 import std.math.traits : isInfinity, isNaN; 87 88 // Special cases. 89 if (isNaN(x) || isInfinity(x)) 90 return x; 91 92 real y = floorImpl(x); 93 if (y < x) 94 y += 1.0; 95 96 return y; 97 } 98 } 99 100 /// 101 @safe pure nothrow @nogc unittest 102 { 103 import std.math.traits : isNaN; 104 105 assert(ceil(+123.456L) == +124); 106 assert(ceil(-123.456L) == -123); 107 assert(ceil(-1.234L) == -1); 108 assert(ceil(-0.123L) == 0); 109 assert(ceil(0.0L) == 0); 110 assert(ceil(+0.123L) == 1); 111 assert(ceil(+1.234L) == 2); 112 assert(ceil(real.infinity) == real.infinity); 113 assert(isNaN(ceil(real.nan))); 114 assert(isNaN(ceil(real.init))); 115 } 116 117 /// ditto 118 double ceil(double x) @trusted pure nothrow @nogc 119 { 120 import std.math.traits : isInfinity, isNaN; 121 122 // Special cases. 123 if (isNaN(x) || isInfinity(x)) 124 return x; 125 126 double y = floorImpl(x); 127 if (y < x) 128 y += 1.0; 129 130 return y; 131 } 132 133 @safe pure nothrow @nogc unittest 134 { 135 import std.math.traits : isNaN; 136 137 assert(ceil(+123.456) == +124); 138 assert(ceil(-123.456) == -123); 139 assert(ceil(-1.234) == -1); 140 assert(ceil(-0.123) == 0); 141 assert(ceil(0.0) == 0); 142 assert(ceil(+0.123) == 1); 143 assert(ceil(+1.234) == 2); 144 assert(ceil(double.infinity) == double.infinity); 145 assert(isNaN(ceil(double.nan))); 146 assert(isNaN(ceil(double.init))); 147 } 148 149 /// ditto 150 float ceil(float x) @trusted pure nothrow @nogc 151 { 152 import std.math.traits : isInfinity, isNaN; 153 154 // Special cases. 155 if (isNaN(x) || isInfinity(x)) 156 return x; 157 158 float y = floorImpl(x); 159 if (y < x) 160 y += 1.0; 161 162 return y; 163 } 164 165 @safe pure nothrow @nogc unittest 166 { 167 import std.math.traits : isNaN; 168 169 assert(ceil(+123.456f) == +124); 170 assert(ceil(-123.456f) == -123); 171 assert(ceil(-1.234f) == -1); 172 assert(ceil(-0.123f) == 0); 173 assert(ceil(0.0f) == 0); 174 assert(ceil(+0.123f) == 1); 175 assert(ceil(+1.234f) == 2); 176 assert(ceil(float.infinity) == float.infinity); 177 assert(isNaN(ceil(float.nan))); 178 assert(isNaN(ceil(float.init))); 179 } 180 181 /************************************** 182 * Returns the value of x rounded downward to the next integer 183 * (toward negative infinity). 184 */ 185 real floor(real x) @trusted pure nothrow @nogc 186 { 187 version (InlineAsm_X87_MSVC) 188 { 189 version (X86_64) 190 { 191 asm pure nothrow @nogc 192 { 193 naked ; 194 fld real ptr [RCX] ; 195 fstcw 8[RSP] ; 196 mov AL,9[RSP] ; 197 mov DL,AL ; 198 and AL,0xC3 ; 199 or AL,0x04 ; // round to -infinity 200 mov 9[RSP],AL ; 201 fldcw 8[RSP] ; 202 frndint ; 203 mov 9[RSP],DL ; 204 fldcw 8[RSP] ; 205 ret ; 206 } 207 } 208 else 209 { 210 short cw; 211 asm pure nothrow @nogc 212 { 213 fld x ; 214 fstcw cw ; 215 mov AL,byte ptr cw+1 ; 216 mov DL,AL ; 217 and AL,0xC3 ; 218 or AL,0x04 ; // round to -infinity 219 mov byte ptr cw+1,AL ; 220 fldcw cw ; 221 frndint ; 222 mov byte ptr cw+1,DL ; 223 fldcw cw ; 224 } 225 } 226 } 227 else 228 { 229 import std.math.traits : isInfinity, isNaN; 230 231 // Special cases. 232 if (isNaN(x) || isInfinity(x) || x == 0.0) 233 return x; 234 235 return floorImpl(x); 236 } 237 } 238 239 /// 240 @safe pure nothrow @nogc unittest 241 { 242 import std.math.traits : isNaN; 243 244 assert(floor(+123.456L) == +123); 245 assert(floor(-123.456L) == -124); 246 assert(floor(+123.0L) == +123); 247 assert(floor(-124.0L) == -124); 248 assert(floor(-1.234L) == -2); 249 assert(floor(-0.123L) == -1); 250 assert(floor(0.0L) == 0); 251 assert(floor(+0.123L) == 0); 252 assert(floor(+1.234L) == 1); 253 assert(floor(real.infinity) == real.infinity); 254 assert(isNaN(floor(real.nan))); 255 assert(isNaN(floor(real.init))); 256 } 257 258 /// ditto 259 double floor(double x) @trusted pure nothrow @nogc 260 { 261 import std.math.traits : isInfinity, isNaN; 262 263 // Special cases. 264 if (isNaN(x) || isInfinity(x) || x == 0.0) 265 return x; 266 267 return floorImpl(x); 268 } 269 270 @safe pure nothrow @nogc unittest 271 { 272 import std.math.traits : isNaN; 273 274 assert(floor(+123.456) == +123); 275 assert(floor(-123.456) == -124); 276 assert(floor(+123.0) == +123); 277 assert(floor(-124.0) == -124); 278 assert(floor(-1.234) == -2); 279 assert(floor(-0.123) == -1); 280 assert(floor(0.0) == 0); 281 assert(floor(+0.123) == 0); 282 assert(floor(+1.234) == 1); 283 assert(floor(double.infinity) == double.infinity); 284 assert(isNaN(floor(double.nan))); 285 assert(isNaN(floor(double.init))); 286 } 287 288 /// ditto 289 float floor(float x) @trusted pure nothrow @nogc 290 { 291 import std.math.traits : isInfinity, isNaN; 292 293 // Special cases. 294 if (isNaN(x) || isInfinity(x) || x == 0.0) 295 return x; 296 297 return floorImpl(x); 298 } 299 300 @safe pure nothrow @nogc unittest 301 { 302 import std.math.traits : isNaN; 303 304 assert(floor(+123.456f) == +123); 305 assert(floor(-123.456f) == -124); 306 assert(floor(+123.0f) == +123); 307 assert(floor(-124.0f) == -124); 308 assert(floor(-1.234f) == -2); 309 assert(floor(-0.123f) == -1); 310 assert(floor(0.0f) == 0); 311 assert(floor(+0.123f) == 0); 312 assert(floor(+1.234f) == 1); 313 assert(floor(float.infinity) == float.infinity); 314 assert(isNaN(floor(float.nan))); 315 assert(isNaN(floor(float.init))); 316 } 317 318 // https://issues.dlang.org/show_bug.cgi?id=6381 319 // floor/ceil should be usable in pure function. 320 @safe pure nothrow unittest 321 { 322 auto x = floor(1.2); 323 auto y = ceil(1.2); 324 } 325 326 /** 327 * Round `val` to a multiple of `unit`. `rfunc` specifies the rounding 328 * function to use; by default this is `rint`, which uses the current 329 * rounding mode. 330 */ 331 Unqual!F quantize(alias rfunc = rint, F)(const F val, const F unit) 332 if (is(typeof(rfunc(F.init)) : F) && isFloatingPoint!F) 333 { 334 import std.math.traits : isInfinity; 335 336 typeof(return) ret = val; 337 if (unit != 0) 338 { 339 const scaled = val / unit; 340 if (!scaled.isInfinity) 341 ret = rfunc(scaled) * unit; 342 } 343 return ret; 344 } 345 346 /// 347 @safe pure nothrow @nogc unittest 348 { 349 import std.math.operations : isClose; 350 351 assert(isClose(12345.6789L.quantize(0.01L), 12345.68L)); 352 assert(isClose(12345.6789L.quantize!floor(0.01L), 12345.67L)); 353 assert(isClose(12345.6789L.quantize(22.0L), 12342.0L)); 354 } 355 356 /// 357 @safe pure nothrow @nogc unittest 358 { 359 import std.math.operations : isClose; 360 import std.math.traits : isNaN; 361 362 assert(isClose(12345.6789L.quantize(0), 12345.6789L)); 363 assert(12345.6789L.quantize(real.infinity).isNaN); 364 assert(12345.6789L.quantize(real.nan).isNaN); 365 assert(real.infinity.quantize(0.01L) == real.infinity); 366 assert(real.infinity.quantize(real.nan).isNaN); 367 assert(real.nan.quantize(0.01L).isNaN); 368 assert(real.nan.quantize(real.infinity).isNaN); 369 assert(real.nan.quantize(real.nan).isNaN); 370 } 371 372 /** 373 * Round `val` to a multiple of `pow(base, exp)`. `rfunc` specifies the 374 * rounding function to use; by default this is `rint`, which uses the 375 * current rounding mode. 376 */ 377 Unqual!F quantize(real base, alias rfunc = rint, F, E)(const F val, const E exp) 378 if (is(typeof(rfunc(F.init)) : F) && isFloatingPoint!F && isIntegral!E) 379 { 380 import std.math.exponential : pow; 381 382 // TODO: Compile-time optimization for power-of-two bases? 383 return quantize!rfunc(val, pow(cast(F) base, exp)); 384 } 385 386 /// ditto 387 Unqual!F quantize(real base, long exp = 1, alias rfunc = rint, F)(const F val) 388 if (is(typeof(rfunc(F.init)) : F) && isFloatingPoint!F) 389 { 390 import std.math.exponential : pow; 391 392 enum unit = cast(F) pow(base, exp); 393 return quantize!rfunc(val, unit); 394 } 395 396 /// 397 @safe pure nothrow @nogc unittest 398 { 399 import std.math.operations : isClose; 400 401 assert(isClose(12345.6789L.quantize!10(-2), 12345.68L)); 402 assert(isClose(12345.6789L.quantize!(10, -2), 12345.68L)); 403 assert(isClose(12345.6789L.quantize!(10, floor)(-2), 12345.67L)); 404 assert(isClose(12345.6789L.quantize!(10, -2, floor), 12345.67L)); 405 406 assert(isClose(12345.6789L.quantize!22(1), 12342.0L)); 407 assert(isClose(12345.6789L.quantize!22, 12342.0L)); 408 } 409 410 @safe pure nothrow @nogc unittest 411 { 412 import std.math.exponential : log10, pow; 413 import std.math.operations : isClose; 414 import std.meta : AliasSeq; 415 416 static foreach (F; AliasSeq!(real, double, float)) 417 {{ 418 const maxL10 = cast(int) F.max.log10.floor; 419 const maxR10 = pow(cast(F) 10, maxL10); 420 assert(isClose((cast(F) 0.9L * maxR10).quantize!10(maxL10), maxR10)); 421 assert(isClose((cast(F)-0.9L * maxR10).quantize!10(maxL10), -maxR10)); 422 423 assert(F.max.quantize(F.min_normal) == F.max); 424 assert((-F.max).quantize(F.min_normal) == -F.max); 425 assert(F.min_normal.quantize(F.max) == 0); 426 assert((-F.min_normal).quantize(F.max) == 0); 427 assert(F.min_normal.quantize(F.min_normal) == F.min_normal); 428 assert((-F.min_normal).quantize(F.min_normal) == -F.min_normal); 429 }} 430 } 431 432 /****************************************** 433 * Rounds x to the nearest integer value, using the current rounding 434 * mode. 435 * 436 * Unlike the rint functions, nearbyint does not raise the 437 * FE_INEXACT exception. 438 */ 439 pragma(inline, true) 440 real nearbyint(real x) @safe pure nothrow @nogc 441 { 442 return core.stdc.math.nearbyintl(x); 443 } 444 445 /// 446 @safe pure unittest 447 { 448 import std.math.traits : isNaN; 449 450 assert(nearbyint(0.4) == 0); 451 assert(nearbyint(0.5) == 0); 452 assert(nearbyint(0.6) == 1); 453 assert(nearbyint(100.0) == 100); 454 455 assert(isNaN(nearbyint(real.nan))); 456 assert(nearbyint(real.infinity) == real.infinity); 457 assert(nearbyint(-real.infinity) == -real.infinity); 458 } 459 460 /********************************** 461 * Rounds x to the nearest integer value, using the current rounding 462 * mode. 463 * 464 * If the return value is not equal to x, the FE_INEXACT 465 * exception is raised. 466 * 467 * $(LREF nearbyint) performs the same operation, but does 468 * not set the FE_INEXACT exception. 469 */ 470 pragma(inline, true) 471 real rint(real x) @safe pure nothrow @nogc 472 { 473 return core.math.rint(x); 474 } 475 ///ditto 476 pragma(inline, true) 477 double rint(double x) @safe pure nothrow @nogc 478 { 479 return core.math.rint(x); 480 } 481 ///ditto 482 pragma(inline, true) 483 float rint(float x) @safe pure nothrow @nogc 484 { 485 return core.math.rint(x); 486 } 487 488 /// 489 @safe unittest 490 { 491 import std.math.traits : isNaN; 492 493 version (IeeeFlagsSupport) resetIeeeFlags(); 494 assert(rint(0.4) == 0); 495 version (IeeeFlagsSupport) assert(ieeeFlags.inexact); 496 497 assert(rint(0.5) == 0); 498 assert(rint(0.6) == 1); 499 assert(rint(100.0) == 100); 500 501 assert(isNaN(rint(real.nan))); 502 assert(rint(real.infinity) == real.infinity); 503 assert(rint(-real.infinity) == -real.infinity); 504 } 505 506 @safe unittest 507 { 508 real function(real) print = &rint; 509 assert(print != null); 510 } 511 512 /*************************************** 513 * Rounds x to the nearest integer value, using the current rounding 514 * mode. 515 * 516 * This is generally the fastest method to convert a floating-point number 517 * to an integer. Note that the results from this function 518 * depend on the rounding mode, if the fractional part of x is exactly 0.5. 519 * If using the default rounding mode (ties round to even integers) 520 * lrint(4.5) == 4, lrint(5.5)==6. 521 */ 522 long lrint(real x) @trusted pure nothrow @nogc 523 { 524 version (InlineAsm_X87) 525 { 526 version (Win64) 527 { 528 asm pure nothrow @nogc 529 { 530 naked; 531 fld real ptr [RCX]; 532 fistp qword ptr 8[RSP]; 533 mov RAX,8[RSP]; 534 ret; 535 } 536 } 537 else 538 { 539 long n; 540 asm pure nothrow @nogc 541 { 542 fld x; 543 fistp n; 544 } 545 return n; 546 } 547 } 548 else 549 { 550 import std.math.traits : floatTraits, RealFormat, MANTISSA_MSB, MANTISSA_LSB; 551 552 alias F = floatTraits!(real); 553 static if (F.realFormat == RealFormat.ieeeDouble) 554 { 555 long result; 556 557 // Rounding limit when casting from real(double) to ulong. 558 enum real OF = 4.50359962737049600000E15L; 559 560 uint* vi = cast(uint*)(&x); 561 562 // Find the exponent and sign 563 uint msb = vi[MANTISSA_MSB]; 564 uint lsb = vi[MANTISSA_LSB]; 565 int exp = ((msb >> 20) & 0x7ff) - 0x3ff; 566 const int sign = msb >> 31; 567 msb &= 0xfffff; 568 msb |= 0x100000; 569 570 if (exp < 63) 571 { 572 if (exp >= 52) 573 result = (cast(long) msb << (exp - 20)) | (lsb << (exp - 52)); 574 else 575 { 576 // Adjust x and check result. 577 const real j = sign ? -OF : OF; 578 x = (j + x) - j; 579 msb = vi[MANTISSA_MSB]; 580 lsb = vi[MANTISSA_LSB]; 581 exp = ((msb >> 20) & 0x7ff) - 0x3ff; 582 msb &= 0xfffff; 583 msb |= 0x100000; 584 585 if (exp < 0) 586 result = 0; 587 else if (exp < 20) 588 result = cast(long) msb >> (20 - exp); 589 else if (exp == 20) 590 result = cast(long) msb; 591 else 592 result = (cast(long) msb << (exp - 20)) | (lsb >> (52 - exp)); 593 } 594 } 595 else 596 { 597 // It is left implementation defined when the number is too large. 598 return cast(long) x; 599 } 600 601 return sign ? -result : result; 602 } 603 else static if (F.realFormat == RealFormat.ieeeExtended || 604 F.realFormat == RealFormat.ieeeExtended53) 605 { 606 long result; 607 608 // Rounding limit when casting from real(80-bit) to ulong. 609 static if (F.realFormat == RealFormat.ieeeExtended) 610 enum real OF = 9.22337203685477580800E18L; 611 else 612 enum real OF = 4.50359962737049600000E15L; 613 614 ushort* vu = cast(ushort*)(&x); 615 uint* vi = cast(uint*)(&x); 616 617 // Find the exponent and sign 618 int exp = (vu[F.EXPPOS_SHORT] & 0x7fff) - 0x3fff; 619 const int sign = (vu[F.EXPPOS_SHORT] >> 15) & 1; 620 621 if (exp < 63) 622 { 623 // Adjust x and check result. 624 const real j = sign ? -OF : OF; 625 x = (j + x) - j; 626 exp = (vu[F.EXPPOS_SHORT] & 0x7fff) - 0x3fff; 627 628 version (LittleEndian) 629 { 630 if (exp < 0) 631 result = 0; 632 else if (exp <= 31) 633 result = vi[1] >> (31 - exp); 634 else 635 result = (cast(long) vi[1] << (exp - 31)) | (vi[0] >> (63 - exp)); 636 } 637 else 638 { 639 if (exp < 0) 640 result = 0; 641 else if (exp <= 31) 642 result = vi[1] >> (31 - exp); 643 else 644 result = (cast(long) vi[1] << (exp - 31)) | (vi[2] >> (63 - exp)); 645 } 646 } 647 else 648 { 649 // It is left implementation defined when the number is too large 650 // to fit in a 64bit long. 651 return cast(long) x; 652 } 653 654 return sign ? -result : result; 655 } 656 else static if (F.realFormat == RealFormat.ieeeQuadruple) 657 { 658 const vu = cast(ushort*)(&x); 659 660 // Find the exponent and sign 661 const sign = (vu[F.EXPPOS_SHORT] >> 15) & 1; 662 if ((vu[F.EXPPOS_SHORT] & F.EXPMASK) - (F.EXPBIAS + 1) > 63) 663 { 664 // The result is left implementation defined when the number is 665 // too large to fit in a 64 bit long. 666 return cast(long) x; 667 } 668 669 // Force rounding of lower bits according to current rounding 670 // mode by adding ±2^-112 and subtracting it again. 671 enum OF = 5.19229685853482762853049632922009600E33L; 672 const j = sign ? -OF : OF; 673 x = (j + x) - j; 674 675 const exp = (vu[F.EXPPOS_SHORT] & F.EXPMASK) - (F.EXPBIAS + 1); 676 const implicitOne = 1UL << 48; 677 auto vl = cast(ulong*)(&x); 678 vl[MANTISSA_MSB] &= implicitOne - 1; 679 vl[MANTISSA_MSB] |= implicitOne; 680 681 long result; 682 683 if (exp < 0) 684 result = 0; 685 else if (exp <= 48) 686 result = vl[MANTISSA_MSB] >> (48 - exp); 687 else 688 result = (vl[MANTISSA_MSB] << (exp - 48)) | (vl[MANTISSA_LSB] >> (112 - exp)); 689 690 return sign ? -result : result; 691 } 692 else 693 { 694 static assert(false, "real type not supported by lrint()"); 695 } 696 } 697 } 698 699 /// 700 @safe pure nothrow @nogc unittest 701 { 702 assert(lrint(4.5) == 4); 703 assert(lrint(5.5) == 6); 704 assert(lrint(-4.5) == -4); 705 assert(lrint(-5.5) == -6); 706 707 assert(lrint(int.max - 0.5) == 2147483646L); 708 assert(lrint(int.max + 0.5) == 2147483648L); 709 assert(lrint(int.min - 0.5) == -2147483648L); 710 assert(lrint(int.min + 0.5) == -2147483648L); 711 } 712 713 static if (real.mant_dig >= long.sizeof * 8) 714 { 715 @safe pure nothrow @nogc unittest 716 { 717 assert(lrint(long.max - 1.5L) == long.max - 1); 718 assert(lrint(long.max - 0.5L) == long.max - 1); 719 assert(lrint(long.min + 0.5L) == long.min); 720 assert(lrint(long.min + 1.5L) == long.min + 2); 721 } 722 } 723 724 /******************************************* 725 * Return the value of x rounded to the nearest integer. 726 * If the fractional part of x is exactly 0.5, the return value is 727 * rounded away from zero. 728 * 729 * Returns: 730 * A `real`. 731 */ 732 auto round(real x) @trusted nothrow @nogc 733 { 734 version (CRuntime_Microsoft) 735 { 736 import std.math.hardware : FloatingPointControl; 737 738 auto old = FloatingPointControl.getControlState(); 739 FloatingPointControl.setControlState( 740 (old & (-1 - FloatingPointControl.roundingMask)) | FloatingPointControl.roundToZero 741 ); 742 x = core.math.rint((x >= 0) ? x + 0.5 : x - 0.5); 743 FloatingPointControl.setControlState(old); 744 return x; 745 } 746 else 747 { 748 return core.stdc.math.roundl(x); 749 } 750 } 751 752 /// 753 @safe nothrow @nogc unittest 754 { 755 assert(round(4.5) == 5); 756 assert(round(5.4) == 5); 757 assert(round(-4.5) == -5); 758 assert(round(-5.1) == -5); 759 } 760 761 // assure purity on Posix 762 version (Posix) 763 { 764 @safe pure nothrow @nogc unittest 765 { 766 assert(round(4.5) == 5); 767 } 768 } 769 770 /********************************************** 771 * Return the value of x rounded to the nearest integer. 772 * 773 * If the fractional part of x is exactly 0.5, the return value is rounded 774 * away from zero. 775 */ 776 long lround(real x) @trusted nothrow @nogc 777 { 778 return core.stdc.math.llroundl(x); 779 } 780 781 /// 782 @safe nothrow @nogc unittest 783 { 784 assert(lround(0.49) == 0); 785 assert(lround(0.5) == 1); 786 assert(lround(1.5) == 2); 787 } 788 789 /** 790 Returns the integer portion of x, dropping the fractional portion. 791 This is also known as "chop" rounding. 792 `pure` on all platforms. 793 */ 794 real trunc(real x) @trusted nothrow @nogc pure 795 { 796 version (InlineAsm_X87_MSVC) 797 { 798 version (X86_64) 799 { 800 asm pure nothrow @nogc 801 { 802 naked ; 803 fld real ptr [RCX] ; 804 fstcw 8[RSP] ; 805 mov AL,9[RSP] ; 806 mov DL,AL ; 807 and AL,0xC3 ; 808 or AL,0x0C ; // round to 0 809 mov 9[RSP],AL ; 810 fldcw 8[RSP] ; 811 frndint ; 812 mov 9[RSP],DL ; 813 fldcw 8[RSP] ; 814 ret ; 815 } 816 } 817 else 818 { 819 short cw; 820 asm pure nothrow @nogc 821 { 822 fld x ; 823 fstcw cw ; 824 mov AL,byte ptr cw+1 ; 825 mov DL,AL ; 826 and AL,0xC3 ; 827 or AL,0x0C ; // round to 0 828 mov byte ptr cw+1,AL ; 829 fldcw cw ; 830 frndint ; 831 mov byte ptr cw+1,DL ; 832 fldcw cw ; 833 } 834 } 835 } 836 else 837 { 838 return core.stdc.math.truncl(x); 839 } 840 } 841 842 /// 843 @safe pure unittest 844 { 845 assert(trunc(0.01) == 0); 846 assert(trunc(0.49) == 0); 847 assert(trunc(0.5) == 0); 848 assert(trunc(1.5) == 1); 849 } 850 851 /***************************************** 852 * Returns x rounded to a long value using the current rounding mode. 853 * If the integer value of x is 854 * greater than long.max, the result is 855 * indeterminate. 856 */ 857 pragma(inline, true) 858 long rndtol(real x) @nogc @safe pure nothrow { return core.math.rndtol(x); } 859 //FIXME 860 ///ditto 861 pragma(inline, true) 862 long rndtol(double x) @safe pure nothrow @nogc { return rndtol(cast(real) x); } 863 //FIXME 864 ///ditto 865 pragma(inline, true) 866 long rndtol(float x) @safe pure nothrow @nogc { return rndtol(cast(real) x); } 867 868 /// 869 @safe unittest 870 { 871 assert(rndtol(1.0) == 1L); 872 assert(rndtol(1.2) == 1L); 873 assert(rndtol(1.7) == 2L); 874 assert(rndtol(1.0001) == 1L); 875 } 876 877 @safe unittest 878 { 879 long function(real) prndtol = &rndtol; 880 assert(prndtol != null); 881 } 882 883 // Helper for floor/ceil 884 T floorImpl(T)(const T x) @trusted pure nothrow @nogc 885 { 886 import std.math.traits : floatTraits, RealFormat; 887 888 alias F = floatTraits!(T); 889 // Take care not to trigger library calls from the compiler, 890 // while ensuring that we don't get defeated by some optimizers. 891 union floatBits 892 { 893 T rv; 894 ushort[T.sizeof/2] vu; 895 896 // Other kinds of extractors for real formats. 897 static if (F.realFormat == RealFormat.ieeeSingle) 898 uint vi; 899 else static if (F.realFormat == RealFormat.ieeeDouble) 900 ulong vi; 901 } 902 floatBits y = void; 903 y.rv = x; 904 905 // Find the exponent (power of 2) 906 // Do this by shifting the raw value so that the exponent lies in the low bits, 907 // then mask out the sign bit, and subtract the bias. 908 static if (F.realFormat == RealFormat.ieeeSingle) 909 { 910 int exp = ((y.vi >> (T.mant_dig - 1)) & 0xff) - 0x7f; 911 enum mantissa_mask = F.MANTISSAMASK_INT; 912 enum sign_shift = 31; 913 } 914 else static if (F.realFormat == RealFormat.ieeeDouble) 915 { 916 long exp = ((y.vi >> (T.mant_dig - 1)) & 0x7ff) - 0x3ff; 917 enum mantissa_mask = F.MANTISSAMASK_LONG; 918 enum sign_shift = 63; 919 } 920 else static if (F.realFormat == RealFormat.ieeeExtended || 921 F.realFormat == RealFormat.ieeeExtended53) 922 { 923 int exp = (y.vu[F.EXPPOS_SHORT] & 0x7fff) - 0x3fff; 924 925 version (LittleEndian) 926 int pos = 0; 927 else 928 int pos = 4; 929 } 930 else static if (F.realFormat == RealFormat.ieeeQuadruple) 931 { 932 int exp = (y.vu[F.EXPPOS_SHORT] & 0x7fff) - 0x3fff; 933 934 version (LittleEndian) 935 int pos = 0; 936 else 937 int pos = 7; 938 } 939 else 940 static assert(false, "Not implemented for this architecture"); 941 942 if (exp < 0) 943 { 944 if (x < 0.0) 945 return -1.0; 946 else 947 return 0.0; 948 } 949 950 static if (F.realFormat == RealFormat.ieeeSingle || 951 F.realFormat == RealFormat.ieeeDouble) 952 { 953 if (exp < (T.mant_dig - 1)) 954 { 955 // Clear all bits representing the fraction part. 956 // Note: the fraction mask represents the floating point number 0.999999... 957 // i.e: `2.0 ^^ (exp - T.mant_dig + 1) * (fraction_mask + 1) == 1.0` 958 const fraction_mask = mantissa_mask >> exp; 959 960 if ((y.vi & fraction_mask) != 0) 961 { 962 // If 'x' is negative, then first substract (1.0 - T.epsilon) from the value. 963 if (y.vi >> sign_shift) 964 y.vi += fraction_mask; 965 y.vi &= ~fraction_mask; 966 } 967 } 968 } 969 else 970 { 971 static if (F.realFormat == RealFormat.ieeeExtended53) 972 exp = (T.mant_dig + 11 - 1) - exp; // mant_dig is really 64 973 else 974 exp = (T.mant_dig - 1) - exp; 975 976 // Zero 16 bits at a time. 977 while (exp >= 16) 978 { 979 version (LittleEndian) 980 y.vu[pos++] = 0; 981 else 982 y.vu[pos--] = 0; 983 exp -= 16; 984 } 985 986 // Clear the remaining bits. 987 if (exp > 0) 988 y.vu[pos] &= 0xffff ^ ((1 << exp) - 1); 989 990 if ((x < 0.0) && (x != y.rv)) 991 y.rv -= 1.0; 992 } 993 994 return y.rv; 995 }