#include <math.h>
#include <stdint.h>
#ifdef ARMLIBM_SET_FLAGS
#include <fenv.h>
#include "required_arithmetic.h"
#pragma STDC FENV_ACCESS ON
double nextafter( double x, double y )
{
union{ double d; uint64_t u;} ux = { x };
uint64_t step = 1;
if( y != y || x != x)
return x + y;
if( y <= x ) {
if( y == x )
return y;
step = -1LL;
}
int64_t signMask = (int64_t) ux.u >> 63;
step = (step ^ signMask) - signMask;
uint64_t absux = ux.u & 0x7fffffffffffffffULL;
if( absux - 0x0010000000000001ULL >= 0x7fefffffffffffffULL - 0x0010000000000001ULL )
{
if( absux == 0ULL ) {
ux.d = y;
required_multiply_double( 0x1.0p-1000, 0x1.0p-1000 );
ux.u = (ux.u & 0x8000000000000000ULL) + 1U;
return ux.d;
}
else if( absux == 0x7fefffffffffffffULL ) {
ux.u += step;
if( 0 == (ux.u & 2ULL) )
{
required_add_double( x, 1.0 ); required_multiply_double( x, x ); }
return ux.d;
}
ux.u += step;
if( 0ULL == (ux.u & 0x7ff0000000000000ULL))
required_multiply_double( 0x1.0p-1000, 0x1.0p-1000 );
return ux.d;
}
ux.u += step;
return ux.d;
}
#else
double nextafter( double x, double y )
{
union{ double d; uint64_t u;} ux = { x };
uint64_t step = 1;
if( y != y || x != x)
return x + y;
if( y <= x ) {
if( y == x )
return y;
step = -1LL;
}
int64_t signMask = (int64_t) ux.u >> 63;
step = (step ^ signMask) - signMask;
uint64_t absux = ux.u & 0x7fffffffffffffffULL;
if( absux == 0ULL ) {
ux.d = y;
ux.u = (ux.u & 0x8000000000000000ULL) + 1U;
return ux.d;
}
ux.u += step;
return ux.d;
}
#endif // ARMLIBM_SET_FLAGS