From owner-freebsd-numerics@FreeBSD.ORG Sat Sep 15 18:21:35 2012 Return-Path: Delivered-To: freebsd-numerics@freebsd.org Received: from mx1.freebsd.org (mx1.freebsd.org [69.147.83.52]) by hub.freebsd.org (Postfix) with ESMTP id 4690D106564A for ; Sat, 15 Sep 2012 18:21:35 +0000 (UTC) (envelope-from stephen@missouri.edu) Received: from wilberforce.math.missouri.edu (wilberforce.math.missouri.edu [128.206.184.213]) by mx1.freebsd.org (Postfix) with ESMTP id DFB7C8FC12 for ; Sat, 15 Sep 2012 18:21:34 +0000 (UTC) Received: from [127.0.0.1] (wilberforce.math.missouri.edu [128.206.184.213]) by wilberforce.math.missouri.edu (8.14.5/8.14.5) with ESMTP id q8FILWjd036494; Sat, 15 Sep 2012 13:21:33 -0500 (CDT) (envelope-from stephen@missouri.edu) Message-ID: <5054C72C.1070903@missouri.edu> Date: Sat, 15 Sep 2012 13:21:32 -0500 From: Stephen Montgomery-Smith User-Agent: Mozilla/5.0 (X11; Linux i686; rv:15.0) Gecko/20120827 Thunderbird/15.0 MIME-Version: 1.0 To: Bruce Evans References: <5048D00B.8010401@missouri.edu> <504D3CCD.2050006@missouri.edu> <504FF726.9060001@missouri.edu> <20120912191556.F1078@besplex.bde.org> <20120912225847.J1771@besplex.bde.org> <50511B40.3070009@missouri.edu> <20120913204808.T1964@besplex.bde.org> <5052A923.9030806@missouri.edu> <20120914143553.X870@besplex.bde.org> <5053CB4F.3030709@missouri.edu> <20120915003408.GA70269@troutmask.apl.washington.edu> <20120915210927.Q2023@besplex.bde.org> In-Reply-To: <20120915210927.Q2023@besplex.bde.org> Content-Type: text/plain; charset=ISO-8859-1; format=flowed Content-Transfer-Encoding: 7bit Cc: freebsd-numerics@freebsd.org, Steve Kargl Subject: Re: cexp error X-BeenThere: freebsd-numerics@freebsd.org X-Mailman-Version: 2.1.5 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, 15 Sep 2012 18:21:35 -0000 On 09/15/2012 07:00 AM, Bruce Evans wrote: > On Fri, 14 Sep 2012, Steve Kargl wrote: > >> On Fri, Sep 14, 2012 at 07:26:55PM -0500, Stephen Montgomery-Smith wrote: >>> On 09/14/2012 01:05 AM, Bruce Evans wrote: >>> >>>> My tests also determined the exact minimum for all multiples up to the >>>> thresholds for "large" multiples in e_rem_pio2.c, except for ld128, >>>> to verify that there are enough bits in the special approximations for >>>> Pi/2 there, except for ld128. The maximum multiples handled there are >>>> 2**28*Pi/2, except for ld128 they are 2**45*Pi/2. 2**45 is too many >>>> to check exhaustively. >>> >>> But presumably one could prove a result to this effect using continued >>> fractions, right? > > Yes. I'm not sure if the theory can produce a good bound for a relatively > small range of x though. > >> See attached. It shows a method for determining the number >> of needed bits. > > It doesn't give the details for continued fractions, but supplies the bound > for 53 mantissa bits (2**-62 -> 61 extra bits), and shows that this was > known in 1992, and describes the 1992 fdlibm better/differently than I > described the current implementation. The 73 (75?) extra that I remembered > must be for 64 mantissa bits. > > I have Kahan's program which does searches using the continued fraction > method (das sent it to me in 2008), but don't really understand it, and > in particular don't know how to hack it to give the bound for the 2**45 > range. Saved results show that I didn't finish checking with it in 2008 > either. I am just making a guess here. But I am thinking that Kahan's program might be based around the following type of argument. First, find integers p, q and r so that |pi - p/q| < 1/r The idea is that r should be much bigger than q. There is a theorem that says you can find a p, q and r so that r=q^2. In any case, I believe continued fractions is the way to find these kinds of numbers. Now suppose that n = the number of bits in the mantissa. Since non-zero multiples of pi/2 are bigger than 1, we are asking the question: how small can |pi - f/2^n| be for any integer f. How use the inequalities: |pi - f/2^n| >= |p/q - f/2^n| - 1/r = |p*2^n - f*q|/(q*2^n) - 1/r. If you have picked r bigger than q*2^n, which we know is possible, then it is just a case of seeing if |p*2^n - f*q| can equal zero. Since we may assume that p and q have no common factors, this can only happen if f=p*2^m and q=2^(n-m). Probably the argument used by Kahan or das@ is sharper than this.