From owner-freebsd-numerics@FreeBSD.ORG Sat Sep 13 06:00:33 2014 Return-Path: Delivered-To: freebsd-numerics@freebsd.org Received: from mx1.freebsd.org (mx1.freebsd.org [8.8.178.115]) (using TLSv1 with cipher ADH-AES256-SHA (256/256 bits)) (No client certificate requested) by hub.freebsd.org (Postfix) with ESMTPS id 08EE68C2 for ; Sat, 13 Sep 2014 06:00:33 +0000 (UTC) Received: from mail106.syd.optusnet.com.au (mail106.syd.optusnet.com.au [211.29.132.42]) by mx1.freebsd.org (Postfix) with ESMTP id BBC26EEE for ; Sat, 13 Sep 2014 06:00:32 +0000 (UTC) Received: from c122-106-147-133.carlnfd1.nsw.optusnet.com.au (c122-106-147-133.carlnfd1.nsw.optusnet.com.au [122.106.147.133]) by mail106.syd.optusnet.com.au (Postfix) with ESMTPS id 3E2503C3083; Sat, 13 Sep 2014 16:00:23 +1000 (EST) Date: Sat, 13 Sep 2014 16:00:22 +1000 (EST) From: Bruce Evans X-X-Sender: bde@besplex.bde.org To: enh Subject: Re: lgamma_r and lgammaf_r return the wrong sign for -0.f In-Reply-To: Message-ID: <20140913152505.E950@besplex.bde.org> References: <20140912220544.GA96714@troutmask.apl.washington.edu> MIME-Version: 1.0 Content-Type: TEXT/PLAIN; charset=US-ASCII; format=flowed X-Optus-CM-Score: 0 X-Optus-CM-Analysis: v=2.1 cv=AOuw8Gd4 c=1 sm=1 tr=0 a=7NqvjVvQucbO2RlWB8PEog==:117 a=PO7r1zJSAAAA:8 a=srA-gRaI36IA:10 a=lreBrSRemiwA:10 a=kj9zAlcOel0A:10 a=JzwRw_2MAAAA:8 a=UaAATAUXzxPjzBn7DBQA:9 a=CjuIK1q_8ugA:10 Cc: freebsd-numerics@freebsd.org, Steve Kargl X-BeenThere: freebsd-numerics@freebsd.org X-Mailman-Version: 2.1.18-1 Precedence: list List-Id: "Discussions of high quality implementation of libm functions." List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , X-List-Received-Date: Sat, 13 Sep 2014 06:00:33 -0000 On Fri, 12 Sep 2014, enh via freebsd-numerics wrote: > On Fri, Sep 12, 2014 at 3:05 PM, Steve Kargl > wrote: >> On Fri, Sep 12, 2014 at 02:15:37PM -0700, enh via freebsd-numerics wrote: >>> if I pass -0.f to lgammaf_r, the sign returned in *signgamp is 1. this >>> is incorrect --- it should be -1. >>> >>> both lgamma_r and lgammaf_r are affected, but the other special cases >>> in those functions look fine to me. >>> >>> this is fixed in OpenBSD and glibc, but FreeBSD and NetBSD both have >>> the same bug. >> >> Are you sure FreeBSD has a bug? From n1256.pdf (committee draft of C99) C99 doesn't even have signgam. But it has the usual complications for +-0 when the sign of the result can be determined for a real function, so it has tgamma(+-0) = +-Inf. POSIX specifies that signgam has the same sign as gamma(), where gamma is spelled with a greek letter. (Bug; it should say the same sign as tgamma(). The one spelled with the greek letter should mean the mathematical one, and that has a complex pole at 0. Complex poles don't have signs.) I looked at what tgamma() does. It is correct for +-0, but for large finite args and +Inf it returns 1/zero. This is only correct for +Inf. >>> patch below (whitespace mangled courtesy of gmail). i'd prefer to wait >>> for this to be fixed in FreeBSD and pull down the fix rather than just >>> fix it locally. The patch can possibly be improved a little by duplicating the zero check ((ix|lx) == 0) instead of duplicating the sign check (hx < 0). This gives the same number of checks, but 1 fewer setting of signgam. In logl(), the checks are arranged more subtly to try to optimize for the usual case. This works best in ld80 long double precision. The single check (hx <= 0) classifies all negative and 0 x. This saves as much as 1 whole cycle or 1.5% on modern x86. It works best in for ld80 since hx doesn't contain any mantissa bits then. When hx contains mantissa bits, its sign bit is usually checked using the C '<' operator as in (hx < 0) and not as a bit. Sometimes the compiler can combine this with another operator and do it as a bit test. The ((ix|lx) == 0) check may be too special for this to be possible except in 64-bit mode where hx, ix, and lx can be treated as a single 64-bit value (then (hlx = 0) works, where int64_t hlx is hx and lx combined). The optimization is going the other way an converting several bit tests into a single signed comparison. Bruce