X-Git-Url: http://nsz.repo.hu/git/?a=blobdiff_plain;f=src%2Fmath%2Fremainder.c;h=6cd089c467cee72945007f92814f60926c98c836;hb=0d0c2f40344640a2a6942dda156509593f51db5d;hp=db176c882164ee5ea790c131c3e8a0fe1d5198f0;hpb=b69f695acedd4ce2798ef9ea28d834ceccc789bd;p=musl diff --git a/src/math/remainder.c b/src/math/remainder.c index db176c88..6cd089c4 100644 --- a/src/math/remainder.c +++ b/src/math/remainder.c @@ -1,70 +1,10 @@ -/* origin: FreeBSD /usr/src/lib/msun/src/e_remainder.c */ -/* - * ==================================================== - * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. - * - * Developed at SunSoft, a Sun Microsystems, Inc. business. - * Permission to use, copy, modify, and distribute this - * software is freely granted, provided that this notice - * is preserved. - * ==================================================== - */ -/* remainder(x,p) - * Return : - * returns x REM p = x - [x/p]*p as if in infinite - * precise arithmetic, where [x/p] is the (infinite bit) - * integer nearest x/p (in half way case choose the even one). - * Method : - * Based on fmod() return x-[x/p]chopped*p exactlp. - */ +#include +#include "libc.h" -#include "libm.h" - -static const double zero = 0.0; - -double remainder(double x, double p) +double remainder(double x, double y) { - int32_t hx,hp; - uint32_t sx,lx,lp; - double p_half; - - EXTRACT_WORDS(hx, lx, x); - EXTRACT_WORDS(hp, lp, p); - sx = hx & 0x80000000; - hp &= 0x7fffffff; - hx &= 0x7fffffff; - - /* purge off exception values */ - if ((hp|lp) == 0) /* p = 0 */ - return (x*p)/(x*p); - if (hx >= 0x7ff00000 || /* x not finite */ - (hp >= 0x7ff00000 && (hp-0x7ff00000 | lp) != 0)) /* p is NaN */ - // FIXME: why long double? - return ((long double)x*p)/((long double)x*p); - - if (hp <= 0x7fdfffff) - x = fmod(x, p+p); /* now x < 2p */ - if (((hx-hp)|(lx-lp)) == 0) - return zero*x; - x = fabs(x); - p = fabs(p); - if (hp < 0x00200000) { - if (x + x > p) { - x -= p; - if (x + x >= p) - x -= p; - } - } else { - p_half = 0.5*p; - if (x > p_half) { - x -= p; - if (x >= p_half) - x -= p; - } - } - GET_HIGH_WORD(hx, x); - if ((hx&0x7fffffff) == 0) - hx = 0; - SET_HIGH_WORD(x, hx^sx); - return x; + int q; + return remquo(x, y, &q); } + +weak_alias(remainder, drem);