system/glibc
3
0
Fork 0
Code Issues Projects Releases Wiki Activity

Optimized remquo for 64-bit platforms

Browse Source
This commit is contained in:
Ulrich Drepper 2011-10-23 21:31:21 -04:00
parent 0275fff8a1
commit bb3129bd25
2 changed files with 114 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
ChangeLog
View File

@ -1,5 +1,7 @@
2011-10-23 Ulrich Drepper <drepper@gmail.com>
* sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c: New file.
* sysdeps/ieee754/dbl-64/e_fmod.c (__ieee754_fmod): Add some branch
prediction.
* sysdeps/ieee754/dbl-64/wordsize-64/e_fmod.c: New file.

112
sysdeps/ieee754/dbl-64/wordsize-64/s_remquo.c Normal file
View File

@ -0,0 +1,112 @@
/* Compute remainder and a congruent to the quotient.
Copyright (C) 1997, 2011 Free Software Foundation, Inc.
This file is part of the GNU C Library.
Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
The GNU C Library is free software; you can redistribute it and/or
modify it under the terms of the GNU Lesser General Public
License as published by the Free Software Foundation; either
version 2.1 of the License, or (at your option) any later version.
The GNU C Library is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
Lesser General Public License for more details.
You should have received a copy of the GNU Lesser General Public
License along with the GNU C Library; if not, write to the Free
Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
02111-1307 USA. */
#include <math.h>
#include "math_private.h"
static const double zero = 0.0;
double
__remquo (double x, double y, int *quo)
{
int64_t hx, hy;
uint64_t sx, qs;
int cquo;
EXTRACT_WORDS64 (hx, x);
EXTRACT_WORDS64 (hy, y);
sx = hx & UINT64_C(0x8000000000000000);
qs = sx ^ (hy & UINT64_C(0x8000000000000000));
hy &= UINT64_C(0x7fffffffffffffff);
hx &= UINT64_C(0x7fffffffffffffff);
/* Purge off exception values. */
if (__builtin_expect (hy == 0, 0))
return (x * y) / (x * y); /* y = 0 */
if (__builtin_expect (hx >= UINT64_C(0x7ff0000000000000) /* x not finite */
|| hy > UINT64_C(0x7ff0000000000000), 0))/* y is NaN */
return (x * y) / (x * y);
if (hy <= UINT64_C(0x7fbfffffffffffff))
x = __ieee754_fmod (x, 8 * y); /* now x < 8y */
if (__builtin_expect (hx == hy, 0))
{
*quo = qs ? -1 : 1;
return zero * x;
}
INSERT_WORDS64 (x, hx);
INSERT_WORDS64 (y, hy);
cquo = 0;
if (x >= 4 * y)
{
x -= 4 * y;
cquo += 4;
}
if (x >= 2 * y)
{
x -= 2 * y;
cquo += 2;
}
if (hy < UINT64_C(0x0020000000000000))
{
if (x + x > y)
{
x -= y;
++cquo;
if (x + x >= y)
{
x -= y;
++cquo;
}
}
}
else
{
double y_half = 0.5 * y;
if (x > y_half)
{
x -= y;
++cquo;
if (x >= y_half)
{
x -= y;
++cquo;
}
}
}
*quo = qs ? -cquo : cquo;
if (sx)
x = -x;
return x;
}
weak_alias (__remquo, remquo)
#ifdef NO_LONG_DOUBLE
strong_alias (__remquo, __remquol)
weak_alias (__remquo, remquol)
#endif