Skip to content

Commit

Permalink
Improve clog, clog10 handling of values with real or imaginary part 1…
Browse files Browse the repository at this point in the history 
… (bug 13629).
  • Loading branch information
Joseph Myers committed Jul 26, 2012
1 parent 3129cfc commit da865e9
Show file tree
Hide file tree
Showing 10 changed files with 318 additions and 54 deletions.
13 changes: 13 additions & 0 deletions ChangeLog
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,18 @@
2012-07-26 Joseph Myers <joseph@codesourcery.com>

[BZ #13629]
* math/s_clog.c (__clog): Use __log1p or direct log1p calculation
if larger part has absolute value 1.0.
* math/s_clog10.c (__clog10): Likewise.
* math/s_clog10f.c (__clog10f): Likewise.
* math/s_clog10l.c (__clog10l): Likewise.
* math/s_clogf.c (__clogf): Likewise.
* math/s_clogl.c (__clogl): Likewise.
* math/libm-test.inc (clog_test): Add more tests.
(clog10_test): Likewise.
* sysdeps/i386/fpu/libm-test-ulps: Update.
* sysdeps/x86_64/fpu/libm-test-ulps: Likewise.

* elf/tst-auditmod1.c (pltenter): Remove all definitions of macro.
(pltexit): Likewise.
(La_regs): Likewise.
Expand Down
40 changes: 40 additions & 0 deletions math/libm-test.inc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2464,6 +2464,26 @@ clog_test (void)
# endif
#endif

TEST_c_c (clog, 1.0L, 0x1.234566p-10L, 6.172834701221959432440126967147726538097e-7L, 1.111110564353742042376451655136933182201e-3L);
TEST_c_c (clog, -1.0L, 0x1.234566p-20L, 5.886877547844618300918562490463748605537e-13L, 3.141591568520436206990380699322226378452L);
TEST_c_c (clog, 0x1.234566p-30L, 1.0L, 5.614163921211322622623353961365728040115e-19L, 1.570796325735258575254858696548386439740L);
TEST_c_c (clog, -0x1.234566p-40L, -1.0L, 5.354083939753840089583620652120903838944e-25L, -1.570796326795931422008642456283782656359L);
TEST_c_c (clog, 0x1.234566p-50L, 1.0L, 5.106052341226425256332038420428899201070e-31L, 1.570796326794895608681734464330528755366L);
TEST_c_c (clog, 0x1.234566p-60L, 1.0L, 4.869510976053643471080816669875627875933e-37L, 1.570796326794896618244456860363082279319L);
TEST_c_c (clog, 0x1p-62L, 1.0L, 2.350988701644575015937473074444491355582e-38L, 1.570796326794896619014481257142650555297L);
TEST_c_c (clog, 0x1p-63L, 1.0L, 5.877471754111437539843682686111228389059e-39L, 1.570796326794896619122901474391200998698L, UNDERFLOW_EXCEPTION_FLOAT);
TEST_c_c (clog, 0x1p-64L, 1.0L, 1.469367938527859384960920671527807097271e-39L, 1.570796326794896619177111583015476220398L, UNDERFLOW_EXCEPTION_FLOAT);
#ifndef TEST_FLOAT
TEST_c_c (clog, 0x1p-510L, 1.0L, 4.450147717014402766180465434664808128438e-308L, 1.570796326794896619231321691639751442099L, UNDERFLOW_EXCEPTION_LDOUBLE_IBM);
TEST_c_c (clog, 0x1p-511L, 1.0L, 1.112536929253600691545116358666202032110e-308L, 1.570796326794896619231321691639751442099L, UNDERFLOW_EXCEPTION_DOUBLE);
TEST_c_c (clog, 0x1p-512L, 1.0L, 2.781342323134001728862790896665505080274e-309L, 1.570796326794896619231321691639751442099L, UNDERFLOW_EXCEPTION_DOUBLE);
#endif
#if defined TEST_LDOUBLE && LDBL_MAX_EXP >= 16384
TEST_c_c (clog, 0x1p-8190L, 1.0L, 6.724206286224187012525355634643505205196e-4932L, 1.570796326794896619231321691639751442099L);
TEST_c_c (clog, 0x1p-8191L, 1.0L, 1.681051571556046753131338908660876301299e-4932L, 1.570796326794896619231321691639751442099L, UNDERFLOW_EXCEPTION);
TEST_c_c (clog, 0x1p-8192L, 1.0L, 4.202628928890116882828347271652190753248e-4933L, 1.570796326794896619231321691639751442099L, UNDERFLOW_EXCEPTION);
#endif

END (clog, complex);
}

Expand Down Expand Up @@ -2593,6 +2613,26 @@ clog10_test (void)
# endif
#endif

TEST_c_c (clog10, 1.0L, 0x1.234566p-10L, 2.680828048441605163181684680300513080769e-7L, 4.825491868832381486767558728169977751564e-4L);
TEST_c_c (clog10, -1.0L, 0x1.234566p-20L, 2.556638434669064077889576526006849923281e-13L, 1.364375882602207106407956770293808181427L);
TEST_c_c (clog10, 0x1.234566p-30L, 1.0L, 2.438200411482400072282924063740535840474e-19L, 6.821881764607257184291586401763604544928e-1L);
TEST_c_c (clog10, -0x1.234566p-40L, -1.0L, 2.325249110681915353442924915876654139373e-25L, -6.821881769213700828789403802671540158935e-1L);
TEST_c_c (clog10, 0x1.234566p-50L, 1.0L, 2.217530356103816369479108963807448194409e-31L, 6.821881769209202348667823902864283966959e-1L);
TEST_c_c (clog10, 0x1.234566p-60L, 1.0L, 2.114801746467415208319767917450504756866e-37L, 6.821881769209206733143018621078368211515e-1L);
TEST_c_c (clog10, 0x1p-61L, 1.0L, 4.084085680564517578238994467153626207224e-38L, 6.821881769209206735545466044044889962925e-1L);
TEST_c_c (clog10, 0x1p-62L, 1.0L, 1.021021420141129394559748616788406551878e-38L, 6.821881769209206736487192085600834406988e-1L, UNDERFLOW_EXCEPTION_FLOAT);
TEST_c_c (clog10, 0x1p-63L, 1.0L, 2.552553550352823486399371541971016379740e-39L, 6.821881769209206736958055106378806629019e-1L, UNDERFLOW_EXCEPTION_FLOAT);
#ifndef TEST_FLOAT
TEST_c_c (clog10, 0x1p-509L, 1.0L, 7.730698388614835910296270976605350994446e-308L, 6.821881769209206737428918127156778851051e-1L, UNDERFLOW_EXCEPTION_LDOUBLE_IBM);
TEST_c_c (clog10, 0x1p-510L, 1.0L, 1.932674597153708977574067744151337748612e-308L, 6.821881769209206737428918127156778851051e-1L, UNDERFLOW_EXCEPTION_DOUBLE);
TEST_c_c (clog10, 0x1p-511L, 1.0L, 4.831686492884272443935169360378344371529e-309L, 6.821881769209206737428918127156778851051e-1L, UNDERFLOW_EXCEPTION_DOUBLE);
#endif
#if defined TEST_LDOUBLE && LDBL_MAX_EXP >= 16384
TEST_c_c (clog10, 0x1p-8189L, 1.0L, 1.168114274114528946314738738025008370069e-4931L, 6.821881769209206737428918127156778851051e-1L);
TEST_c_c (clog10, 0x1p-8190L, 1.0L, 2.920285685286322365786846845062520925172e-4932L, 6.821881769209206737428918127156778851051e-1L, UNDERFLOW_EXCEPTION);
TEST_c_c (clog10, 0x1p-8191L, 1.0L, 7.300714213215805914467117112656302312931e-4933L, 6.821881769209206737428918127156778851051e-1L, UNDERFLOW_EXCEPTION);
#endif

END (clog10, complex);
}

Expand Down
36 changes: 27 additions & 9 deletions math/s_clog.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -41,31 +41,49 @@ __clog (__complex__ double x)
{
/* Neither real nor imaginary part is NaN. */
double absx = fabs (__real__ x), absy = fabs (__imag__ x);
double d;
int scale = 0;

if (absx < absy)
{
double t = absx;
absx = absy;
absy = t;
}

if (absx > DBL_MAX / 2.0)
{
scale = -1;
absx = __scalbn (absx, scale);
absy = (absy >= DBL_MIN * 2.0 ? __scalbn (absy, scale) : 0.0);
}
else if (absy > DBL_MAX / 2.0)
{
scale = -1;
absx = (absx >= DBL_MIN * 2.0 ? __scalbn (absx, scale) : 0.0);
absy = __scalbn (absy, scale);
}
else if (absx < DBL_MIN && absy < DBL_MIN)
{
scale = DBL_MANT_DIG;
absx = __scalbn (absx, scale);
absy = __scalbn (absy, scale);
}

d = __ieee754_hypot (absx, absy);
if (absx == 1.0 && scale == 0)
{
double absy2 = absy * absy;
if (absy2 <= DBL_MIN * 2.0)
{
#if __FLT_EVAL_METHOD__ == 0
__real__ result = absy2 / 2.0 - absy2 * absy2 / 4.0;
#else
volatile double force_underflow = absy2 * absy2 / 4.0;
__real__ result = absy2 / 2.0 - force_underflow;
#endif
}
else
__real__ result = __log1p (absy2) / 2.0;
}
else
{
double d = __ieee754_hypot (absx, absy);
__real__ result = __ieee754_log (d) - scale * M_LN2;
}

__real__ result = __ieee754_log (d) - scale * M_LN2;
__imag__ result = __ieee754_atan2 (__imag__ x, __real__ x);
}
else
Expand Down
36 changes: 27 additions & 9 deletions math/s_clog10.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,31 +44,49 @@ __clog10 (__complex__ double x)
{
/* Neither real nor imaginary part is NaN. */
double absx = fabs (__real__ x), absy = fabs (__imag__ x);
double d;
int scale = 0;

if (absx < absy)
{
double t = absx;
absx = absy;
absy = t;
}

if (absx > DBL_MAX / 2.0)
{
scale = -1;
absx = __scalbn (absx, scale);
absy = (absy >= DBL_MIN * 2.0 ? __scalbn (absy, scale) : 0.0);
}
else if (absy > DBL_MAX / 2.0)
{
scale = -1;
absx = (absx >= DBL_MIN * 2.0 ? __scalbn (absx, scale) : 0.0);
absy = __scalbn (absy, scale);
}
else if (absx < DBL_MIN && absy < DBL_MIN)
{
scale = DBL_MANT_DIG;
absx = __scalbn (absx, scale);
absy = __scalbn (absy, scale);
}

d = __ieee754_hypot (absx, absy);
if (absx == 1.0 && scale == 0)
{
double absy2 = absy * absy;
if (absy2 <= DBL_MIN * 2.0 * M_LN10)
{
#if __FLT_EVAL_METHOD__ == 0
__real__ result = (absy2 / 2.0 - absy2 * absy2 / 4.0) * M_LOG10E;
#else
volatile double force_underflow = absy2 * absy2 / 4.0;
__real__ result = (absy2 / 2.0 - force_underflow) * M_LOG10E;
#endif
}
else
__real__ result = __log1p (absy2) * (M_LOG10E / 2.0);
}
else
{
double d = __ieee754_hypot (absx, absy);
__real__ result = __ieee754_log10 (d) - scale * M_LOG10_2;
}

__real__ result = __ieee754_log10 (d) - scale * M_LOG10_2;
__imag__ result = M_LOG10E * __ieee754_atan2 (__imag__ x, __real__ x);
}
else
Expand Down
38 changes: 29 additions & 9 deletions math/s_clog10f.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,31 +44,51 @@ __clog10f (__complex__ float x)
{
/* Neither real nor imaginary part is NaN. */
float absx = fabsf (__real__ x), absy = fabsf (__imag__ x);
float d;
int scale = 0;

if (absx < absy)
{
float t = absx;
absx = absy;
absy = t;
}

if (absx > FLT_MAX / 2.0f)
{
scale = -1;
absx = __scalbnf (absx, scale);
absy = (absy >= FLT_MIN * 2.0f ? __scalbnf (absy, scale) : 0.0f);
}
else if (absy > FLT_MAX / 2.0f)
{
scale = -1;
absx = (absx >= FLT_MIN * 2.0f ? __scalbnf (absx, scale) : 0.0f);
absy = __scalbnf (absy, scale);
}
else if (absx < FLT_MIN && absy < FLT_MIN)
{
scale = FLT_MANT_DIG;
absx = __scalbnf (absx, scale);
absy = __scalbnf (absy, scale);
}

d = __ieee754_hypotf (absx, absy);
if (absx == 1.0f && scale == 0)
{
float absy2 = absy * absy;
if (absy2 <= FLT_MIN * 2.0f * (float) M_LN10)
{
#if __FLT_EVAL_METHOD__ == 0
__real__ result
= (absy2 / 2.0f - absy2 * absy2 / 4.0f) * (float) M_LOG10E;
#else
volatile float force_underflow = absy2 * absy2 / 4.0f;
__real__ result
= (absy2 / 2.0f - force_underflow) * (float) M_LOG10E;
#endif
}
else
__real__ result = __log1pf (absy2) * ((float) M_LOG10E / 2.0f);
}
else
{
float d = __ieee754_hypotf (absx, absy);
__real__ result = __ieee754_log10f (d) - scale * M_LOG10_2f;
}

__real__ result = __ieee754_log10f (d) - scale * M_LOG10_2f;
__imag__ result = M_LOG10E * __ieee754_atan2f (__imag__ x, __real__ x);
}
else
Expand Down
30 changes: 21 additions & 9 deletions math/s_clog10l.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -44,31 +44,43 @@ __clog10l (__complex__ long double x)
{
/* Neither real nor imaginary part is NaN. */
long double absx = fabsl (__real__ x), absy = fabsl (__imag__ x);
long double d;
int scale = 0;

if (absx < absy)
{
long double t = absx;
absx = absy;
absy = t;
}

if (absx > LDBL_MAX / 2.0L)
{
scale = -1;
absx = __scalbnl (absx, scale);
absy = (absy >= LDBL_MIN * 2.0L ? __scalbnl (absy, scale) : 0.0L);
}
else if (absy > LDBL_MAX / 2.0L)
{
scale = -1;
absx = (absx >= LDBL_MIN * 2.0L ? __scalbnl (absx, scale) : 0.0L);
absy = __scalbnl (absy, scale);
}
else if (absx < LDBL_MIN && absy < LDBL_MIN)
{
scale = LDBL_MANT_DIG;
absx = __scalbnl (absx, scale);
absy = __scalbnl (absy, scale);
}

d = __ieee754_hypotl (absx, absy);
if (absx == 1.0L && scale == 0)
{
long double absy2 = absy * absy;
if (absy2 <= LDBL_MIN * 2.0L * M_LN10l)
__real__ result
= (absy2 / 2.0L - absy2 * absy2 / 4.0L) * M_LOG10El;
else
__real__ result = __log1pl (absy2) * (M_LOG10El / 2.0L);
}
else
{
long double d = __ieee754_hypotl (absx, absy);
__real__ result = __ieee754_log10l (d) - scale * M_LOG10_2l;
}

__real__ result = __ieee754_log10l (d) - scale * M_LOG10_2l;
__imag__ result = M_LOG10El * __ieee754_atan2l (__imag__ x, __real__ x);
}
else
Expand Down
36 changes: 27 additions & 9 deletions math/s_clogf.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -41,31 +41,49 @@ __clogf (__complex__ float x)
{
/* Neither real nor imaginary part is NaN. */
float absx = fabsf (__real__ x), absy = fabsf (__imag__ x);
float d;
int scale = 0;

if (absx < absy)
{
float t = absx;
absx = absy;
absy = t;
}

if (absx > FLT_MAX / 2.0f)
{
scale = -1;
absx = __scalbnf (absx, scale);
absy = (absy >= FLT_MIN * 2.0f ? __scalbnf (absy, scale) : 0.0f);
}
else if (absy > FLT_MAX / 2.0f)
{
scale = -1;
absx = (absx >= FLT_MIN * 2.0f ? __scalbnf (absx, scale) : 0.0f);
absy = __scalbnf (absy, scale);
}
else if (absx < FLT_MIN && absy < FLT_MIN)
{
scale = FLT_MANT_DIG;
absx = __scalbnf (absx, scale);
absy = __scalbnf (absy, scale);
}

d = __ieee754_hypotf (absx, absy);
if (absx == 1.0f && scale == 0)
{
float absy2 = absy * absy;
if (absy2 <= FLT_MIN * 2.0f)
{
#if __FLT_EVAL_METHOD__ == 0
__real__ result = absy2 / 2.0f - absy2 * absy2 / 4.0f;
#else
volatile float force_underflow = absy2 * absy2 / 4.0f;
__real__ result = absy2 / 2.0f - force_underflow;
#endif
}
else
__real__ result = __log1pf (absy2) / 2.0f;
}
else
{
float d = __ieee754_hypotf (absx, absy);
__real__ result = __ieee754_logf (d) - scale * (float) M_LN2;
}

__real__ result = __ieee754_logf (d) - scale * (float) M_LN2;
__imag__ result = __ieee754_atan2f (__imag__ x, __real__ x);
}
else
Expand Down
Loading

0 comments on commit da865e9

Please sign in to comment.