1 /*
   2  * CDDL HEADER START
   3  *
   4  * The contents of this file are subject to the terms of the
   5  * Common Development and Distribution License (the "License").
   6  * You may not use this file except in compliance with the License.
   7  *
   8  * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
   9  * or http://www.opensolaris.org/os/licensing.
  10  * See the License for the specific language governing permissions
  11  * and limitations under the License.
  12  *
  13  * When distributing Covered Code, include this CDDL HEADER in each
  14  * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
  15  * If applicable, add the following below this CDDL HEADER, with the
  16  * fields enclosed by brackets "[]" replaced with your own identifying
  17  * information: Portions Copyright [yyyy] [name of copyright owner]
  18  *
  19  * CDDL HEADER END
  20  */
  21 
  22 /*
  23  * Copyright 2011 Nexenta Systems, Inc.  All rights reserved.
  24  */
  25 
  26 /*
  27  * Copyright 2006 Sun Microsystems, Inc.  All rights reserved.
  28  * Use is subject to license terms.
  29  */
  30 
  31 #pragma weak nexttoward = __nexttoward
  32 
  33 /*
  34  * nexttoward(x, y) delivers the next representable number after x
  35  * in the direction of y.  If x and y are both zero, the result is
  36  * zero with the same sign as y.  If either x or y is NaN, the result
  37  * is NaN.
  38  *
  39  * If x != y and the result is infinite, overflow is raised; if
  40  * x != y and the result is subnormal or zero, underflow is raised.
  41  * (This is wrong, but it's what C99 apparently wants.)
  42  */
  43 
  44 #include "libm.h"
  45 
  46 #if defined(__sparc)
  47 static union {
  48         unsigned i[2];
  49         double d;
  50 } C[] = {
  51         0x00100000, 0, 0x7fe00000, 0, 0x7fffffff, 0xffffffff
  52 };
  53 
  54 #define tiny            C[0].d
  55 #define huge            C[1].d
  56 #define qnan            C[2].d
  57 
  58 enum fcc_type {
  59         fcc_equal = 0, fcc_less = 1, fcc_greater = 2, fcc_unordered = 3
  60 };
  61 
  62 #ifdef __sparcv9
  63 #define _Q_cmp          _Qp_cmp
  64 #endif
  65 
  66 extern enum fcc_type _Q_cmp(const long double *, const long double *);
  67 
  68 double
  69 __nexttoward(double x, long double y)
  70 {
  71         union {
  72                 unsigned i[2];
  73                 double d;
  74         } xx;
  75         union {
  76                 unsigned i[4];
  77                 long double q;
  78         } yy;
  79 
  80         long double lx;
  81         unsigned hx;
  82         volatile double dummy;
  83         enum fcc_type rel;
  84 
  85         /*
  86          * It would be somewhat more efficient to check for NaN and
  87          * zero operands before converting x to long double and then
  88          * to code the comparison in line rather than calling _Q_cmp.
  89          * However, since this code probably won't get used much,
  90          * I'm opting in favor of simplicity instead.
  91          */
  92         lx = xx.d = x;
  93         hx = (xx.i[0] & ~0x80000000) | xx.i[1];
  94 
  95         /* check for each of four possible orderings */
  96         rel = _Q_cmp(&lx, &y);
  97 
  98         if (rel == fcc_unordered)
  99                 return (qnan);
 100 
 101         if (rel == fcc_equal) {
 102                 if (hx == 0) {  /* x is zero; return zero with y's sign */
 103                         yy.q = y;
 104                         xx.i[0] = yy.i[0];
 105                         return (xx.d);
 106                 }
 107 
 108                 return (x);
 109         }
 110 
 111         if (rel == fcc_less) {
 112                 if (hx == 0) {                  /* x is zero */
 113                         xx.i[0] = 0;
 114                         xx.i[1] = 0x00000001;
 115                 } else if ((int)xx.i[0] >= 0) {      /* x is positive */
 116                         if (++xx.i[1] == 0)
 117                                 xx.i[0]++;
 118                 } else {
 119                         if (xx.i[1]-- == 0)
 120                                 xx.i[0]--;
 121                 }
 122         } else {
 123                 if (hx == 0) {                  /* x is zero */
 124                         xx.i[0] = 0x80000000;
 125                         xx.i[1] = 0x00000001;
 126                 } else if ((int)xx.i[0] >= 0) {      /* x is positive */
 127                         if (xx.i[1]-- == 0)
 128                                 xx.i[0]--;
 129                 } else {
 130                         if (++xx.i[1] == 0)
 131                                 xx.i[0]++;
 132                 }
 133         }
 134 
 135         /* raise exceptions as needed */
 136         hx = xx.i[0] & ~0x80000000;
 137 
 138         if (hx == 0x7ff00000) {
 139                 dummy = huge;
 140                 dummy *= huge;
 141         } else if (hx < 0x00100000) {
 142                 dummy = tiny;
 143                 dummy *= tiny;
 144         }
 145 
 146         return (xx.d);
 147 }
 148 #elif defined(__x86)
 149 static union {
 150         unsigned i[2];
 151         double d;
 152 } C[] = {
 153         0, 0x00100000, 0, 0x7fe00000, };
 154 
 155 #define tiny            C[0].d
 156 #define huge            C[1].d
 157 
 158 double
 159 __nexttoward(double x, long double y)
 160 {
 161         union {
 162                 unsigned i[2];
 163                 double d;
 164         } xx;
 165 
 166         unsigned hx;
 167         long double lx;
 168         volatile double dummy;
 169 
 170         lx = xx.d = x;
 171         hx = (xx.i[1] & ~0x80000000) | xx.i[0];
 172 
 173         /* check for each of four possible orderings */
 174         if (isunordered(lx, y))
 175                 return ((double)(lx + y));
 176 
 177         if (lx == y)
 178                 return ((double)y);
 179 
 180         if (lx < y) {
 181                 if (hx == 0) {                  /* x is zero */
 182                         xx.i[0] = 0x00000001;
 183                         xx.i[1] = 0;
 184                 } else if ((int)xx.i[1] >= 0) {      /* x is positive */
 185                         if (++xx.i[0] == 0)
 186                                 xx.i[1]++;
 187                 } else {
 188                         if (xx.i[0]-- == 0)
 189                                 xx.i[1]--;
 190                 }
 191         } else {
 192                 if (hx == 0) {                  /* x is zero */
 193                         xx.i[0] = 0x00000001;
 194                         xx.i[1] = 0x80000000;
 195                 } else if ((int)xx.i[1] >= 0) {      /* x is positive */
 196                         if (xx.i[0]-- == 0)
 197                                 xx.i[1]--;
 198                 } else {
 199                         if (++xx.i[0] == 0)
 200                                 xx.i[1]++;
 201                 }
 202         }
 203 
 204         /* raise exceptions as needed */
 205         hx = xx.i[1] & ~0x80000000;
 206 
 207         if (hx == 0x7ff00000) {
 208                 dummy = huge;
 209                 dummy *= huge;
 210         } else if (hx < 0x00100000) {
 211                 dummy = tiny;
 212                 dummy *= tiny;
 213         }
 214 
 215         return (xx.d);
 216 }
 217 #else
 218 #error Unknown architecture
 219 #endif