Subject: lib/22599: math library i387 log1p algorithm too simplistic
To: None <gnats-bugs@gnats.netbsd.org>
From: None <ray@mcs.vuw.ac.nz>
List: netbsd-bugs
Date: 08/26/2003 13:28:56
>Number:         22599
>Category:       lib
>Synopsis:       i387 assembly version of log1p too simplistic
>Confidential:   no
>Severity:       serious
>Priority:       medium
>Responsible:    lib-bug-people
>State:          open
>Class:          sw-bug
>Submitter-Id:   net
>Arrival-Date:   Tue Aug 26 01:29:00 UTC 2003
>Closed-Date:
>Last-Modified:
>Originator:     Ray Brownrigg
>Release:        NetBSD 1.6W
>Organization:
Victoria University of Wellington
>Environment:
System: NetBSD turakirae.mcs.vuw.ac.nz 1.6W NetBSD 1.6W (MCS_WORKSTATION) #1: Wed Aug 20 17:33:27 NZST 2003 mark@dellc64h2.mcs.vuw.ac.nz:/mnt/SAVE/build.obj/mnt/src/src/sys/arch/i386/compile/MCS_WORKSTATION i386
Architecture: i386
Machine: i386
>Description:
	The i387 version of the log1p math library function erroneously
	shuns the appropriate instruction, and reverts to a simplistic
	algorithm that defeats the whole purpose of the function.
>How-To-Repeat:
turakirae>  cat bug-log1p.c
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

int
main(int argc, char* argv[])
{
    int k;
    double x;

    for (k = 10; k <= 20; ++k)
    {
        x = pow(10, -k);
        printf("%3d\t%g\n", k, log1p(x));
    }

    return (EXIT_SUCCESS);
}
turakirae>  cc bug-log1p.c -lm && ./a.out
 10     1e-10
 11     1e-11
 12     1.00009e-12
 13     9.99201e-14
 14     9.99201e-15
 15     1.11022e-15
 16     0
 17     0
 18     0
 19     0
 20     0
turakirae>  # Using the non-i387 code works fine:
turakirae>  cc -c /src/work/src/lib/libm/src/s_log1p.c
turakirae>  cc bug-log1p.c s_log1p.o -lm && ./a.out
 10     1e-10
 11     1e-11
 12     1e-12
 13     1e-13
 14     1e-14
 15     1e-15
 16     1e-16
 17     1e-17
 18     1e-18
 19     1e-19
 20     1e-20
turakirae>

>Fix:
Either:
1) dispense with the i387 version of log1p and just use
.../src/lib/libm/src/s_log1p.c or:
2) modify .../src/lib/libm/arch/i387/s_log1p.S to properly use the
fyl2xp1 instruction when the argument is in the appropriate but
restricted (in some older CPUs) range.
>Release-Note:
>Audit-Trail:
>Unformatted: