s_nextafterf.c 2.3 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596
  1. /* s_nextafterf.c -- float version of s_nextafter.c.
  2. * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
  3. */
  4. /*
  5. * ====================================================
  6. * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
  7. *
  8. * Developed at SunPro, a Sun Microsystems, Inc. business.
  9. * Permission to use, copy, modify, and distribute this
  10. * software is freely granted, provided that this notice
  11. * is preserved.
  12. * ====================================================
  13. */
  14. #include "math.h"
  15. #include "math_private.h"
  16. float nextafterf(float x, float y)
  17. {
  18. int32_t hx, hy, ix, iy;
  19. GET_FLOAT_WORD(hx, x);
  20. GET_FLOAT_WORD(hy, y);
  21. ix = hx & 0x7fffffff; /* |x| */
  22. iy = hy & 0x7fffffff; /* |y| */
  23. /* x is nan or y is nan? */
  24. if ((ix > 0x7f800000) || (iy > 0x7f800000))
  25. return x + y;
  26. if (x == y)
  27. return y;
  28. if (ix == 0) { /* x == 0? */
  29. /* glibc 2.4 does not seem to set underflow? */
  30. /* float u; */
  31. /* return +-minsubnormal */
  32. SET_FLOAT_WORD(x, (hy & 0x80000000) | 1);
  33. /* u = x * x; raise underflow flag */
  34. /* math_force_eval(u); */
  35. return x;
  36. }
  37. if (hx >= 0) { /* x > 0 */
  38. if (hx > hy) { /* x > y: x -= ulp */
  39. hx -= 1;
  40. } else { /* x < y: x += ulp */
  41. hx += 1;
  42. }
  43. } else { /* x < 0 */
  44. if (hy >= 0 || hx > hy) { /* x < y: x -= ulp */
  45. hx -= 1;
  46. } else { /* x > y: x += ulp */
  47. hx += 1;
  48. }
  49. }
  50. hy = hx & 0x7f800000;
  51. if (hy >= 0x7f800000) {
  52. x = x + x; /* overflow */
  53. return x; /* overflow */
  54. }
  55. if (hy < 0x00800000) {
  56. float u = x * x; /* underflow */
  57. math_force_eval(u); /* raise underflow flag */
  58. }
  59. SET_FLOAT_WORD(x, hx);
  60. return x;
  61. }
  62. #if 0
  63. /* "testprog N a b"
  64. * calculates a = nextafterf(a, b) and prints a as float
  65. * and as raw bytes; repeats it N times.
  66. */
  67. #include <stdio.h>
  68. #include <stdlib.h>
  69. #include <math.h>
  70. int main(int argc, char **argv)
  71. {
  72. int cnt, i;
  73. float a, b;
  74. cnt = atoi(argv[1]);
  75. a = strtod(argv[2], NULL);
  76. b = strtod(argv[3], NULL);
  77. while (cnt-- > 0) {
  78. for (i = 0; i < sizeof(a); i++) {
  79. unsigned char c = ((char*)(&a))[i];
  80. printf("%x%x", (c >> 4), (c & 0xf));
  81. }
  82. printf(" %f\n", a);
  83. a = nextafterf(a, b);
  84. }
  85. return 0;
  86. }
  87. #endif