1  /* @(#)s_nextafter.c 5.1 93/09/24 */
       2  /*
       3   * ====================================================
       4   * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
       5   *
       6   * Developed at SunPro, a Sun Microsystems, Inc. business.
       7   * Permission to use, copy, modify, and distribute this
       8   * software is freely granted, provided that this notice
       9   * is preserved.
      10   * ====================================================
      11   */
      12  
      13  #if defined(LIBM_SCCS) && !defined(lint)
      14  static char rcsid[] = "$NetBSD: s_nextafter.c,v 1.8 1995/05/10 20:47:58 jtc Exp $";
      15  #endif
      16  
      17  /* IEEE functions
      18   *	nextafter(x,y)
      19   *	return the next machine floating-point number of x in the
      20   *	direction toward y.
      21   *   Special cases:
      22   */
      23  
      24  /* Ugly hack so that the aliasing works.  */
      25  #define __nexttoward __internal___nexttoward
      26  #define nexttoward __internal_nexttoward
      27  
      28  #include <errno.h>
      29  #include <math.h>
      30  #include <math-barriers.h>
      31  #include <math_private.h>
      32  #include <float.h>
      33  #include <libm-alias-double.h>
      34  
      35  double __nextafter(double x, double y)
      36  {
      37  	int32_t hx,hy,ix,iy;
      38  	uint32_t lx,ly;
      39  
      40  	EXTRACT_WORDS(hx,lx,x);
      41  	EXTRACT_WORDS(hy,ly,y);
      42  	ix = hx&0x7fffffff;		/* |x| */
      43  	iy = hy&0x7fffffff;		/* |y| */
      44  
      45  	if(((ix>=0x7ff00000)&&((ix-0x7ff00000)|lx)!=0) ||   /* x is nan */
      46  	   ((iy>=0x7ff00000)&&((iy-0x7ff00000)|ly)!=0))     /* y is nan */
      47  	   return x+y;
      48  	if(x==y) return y;		/* x=y, return y */
      49  	if((ix|lx)==0) {			/* x == 0 */
      50  	    double u;
      51  	    INSERT_WORDS(x,hy&0x80000000,1);	/* return +-minsubnormal */
      52  	    u = math_opt_barrier (x);
      53  	    u = u*u;
      54  	    math_force_eval (u);		/* raise underflow flag */
      55  	    return x;
      56  	}
      57  	if(hx>=0) {				/* x > 0 */
      58  	    if(hx>hy||((hx==hy)&&(lx>ly))) {	/* x > y, x -= ulp */
      59  		if(lx==0) hx -= 1;
      60  		lx -= 1;
      61  	    } else {				/* x < y, x += ulp */
      62  		lx += 1;
      63  		if(lx==0) hx += 1;
      64  	    }
      65  	} else {				/* x < 0 */
      66  	    if(hy>=0||hx>hy||((hx==hy)&&(lx>ly))){/* x < y, x -= ulp */
      67  		if(lx==0) hx -= 1;
      68  		lx -= 1;
      69  	    } else {				/* x > y, x += ulp */
      70  		lx += 1;
      71  		if(lx==0) hx += 1;
      72  	    }
      73  	}
      74  	hy = hx&0x7ff00000;
      75  	if(hy>=0x7ff00000) {
      76  	  double u = x+x;	/* overflow  */
      77  	  math_force_eval (u);
      78  	  __set_errno (ERANGE);
      79  	}
      80  	if(hy<0x00100000) {
      81  	    double u = x*x;			/* underflow */
      82  	    math_force_eval (u);		/* raise underflow flag */
      83  	  __set_errno (ERANGE);
      84  	}
      85  	INSERT_WORDS(x,hx,lx);
      86  	return x;
      87  }
      88  libm_alias_double (__nextafter, nextafter)
      89  #ifdef NO_LONG_DOUBLE
      90  strong_alias (__nextafter, __nexttowardl)
      91  weak_alias (__nexttowardl, nexttowardl)
      92  #undef __nexttoward
      93  strong_alias (__nextafter, __nexttoward)
      94  #undef nexttoward
      95  weak_alias (__nextafter, nexttoward)
      96  #endif