From owner-freebsd-standards@FreeBSD.ORG Thu Dec 6 23:15:44 2007 Return-Path: Delivered-To: freebsd-standards@FreeBSD.ORG Received: from mx1.freebsd.org (mx1.freebsd.org [IPv6:2001:4f8:fff6::34]) by hub.freebsd.org (Postfix) with ESMTP id 3AE6E16A419 for ; Thu, 6 Dec 2007 23:15:44 +0000 (UTC) (envelope-from sgk@troutmask.apl.washington.edu) Received: from troutmask.apl.washington.edu (troutmask.apl.washington.edu [128.208.78.105]) by mx1.freebsd.org (Postfix) with ESMTP id 07E8313C442 for ; Thu, 6 Dec 2007 23:15:44 +0000 (UTC) (envelope-from sgk@troutmask.apl.washington.edu) Received: from troutmask.apl.washington.edu (localhost.apl.washington.edu [127.0.0.1]) by troutmask.apl.washington.edu (8.14.1/8.14.1) with ESMTP id lB6NBhpP064497 for ; Thu, 6 Dec 2007 15:11:43 -0800 (PST) (envelope-from sgk@troutmask.apl.washington.edu) Received: (from sgk@localhost) by troutmask.apl.washington.edu (8.14.1/8.14.1/Submit) id lB6NBh9I064496 for freebsd-standards@FreeBSD.ORG; Thu, 6 Dec 2007 15:11:43 -0800 (PST) (envelope-from sgk) Date: Thu, 6 Dec 2007 15:11:43 -0800 From: Steve Kargl To: freebsd-standards@FreeBSD.ORG Message-ID: <20071206231143.GA63969@troutmask.apl.washington.edu> References: <20071012180959.GA36345@troutmask.apl.washington.edu> <20071206090833.GA95428@VARK.MIT.EDU> Mime-Version: 1.0 Content-Type: text/plain; charset=us-ascii Content-Disposition: inline In-Reply-To: <20071206090833.GA95428@VARK.MIT.EDU> User-Agent: Mutt/1.4.2.3i Cc: Subject: Re: [PATCH] hypotl, cabsl, and code removal in cabs X-BeenThere: freebsd-standards@freebsd.org X-Mailman-Version: 2.1.5 Precedence: list List-Id: Standards compliance List-Unsubscribe: , List-Archive: List-Post: List-Help: List-Subscribe: , X-List-Received-Date: Thu, 06 Dec 2007 23:15:44 -0000 On Thu, Dec 06, 2007 at 04:08:33AM -0500, David Schultz wrote: > I like the approach. One concern, though, is that IEEE 754 > requires sqrt to always be correctly rounded, and one might > reasonably expect the same from hypot. If we believe that, then > the intermediate quantity a*a + b*b needs to be computed exactly; > if it is rounded before we take the square root, the double > rounding may cause us to round the wrong way. This basically > requires computing things with twice as many bits of precision. > Arguably, though, your code is correct in the most important > ways (avoiding underflow/overflow). > > Also, umm, I've been busy and unable to pay attention for a while, > so forgive me if I'm missing something, but isn't it the case that > we don't have a sqrtl(), except for the gcc builtin on some > architectures? David, bde pointed me to the right file in src/libm/ieee that explains the rounding issues with hypotl. I haven't had a chance to update my implementation to use extra care in the evaluation of a*a+b*b. As to the sqrtl question, I have an implementation that supposely does correct rounding in all rounding modes. It is restricted to 64-bit significand long doubles. The code does not use bit twiddle; instead, it uses fenv. -- Steve