From mboxrd@z Thu Jan 1 00:00:00 1970 X-Spam-Checker-Version: SpamAssassin 3.4.4 (2020-01-24) on inbox.vuxu.org X-Spam-Level: X-Spam-Status: No, score=-3.3 required=5.0 tests=MAILING_LIST_MULTI, RCVD_IN_DNSWL_MED,RCVD_IN_MSPIKE_H3,RCVD_IN_MSPIKE_WL autolearn=ham autolearn_force=no version=3.4.4 Received: (qmail 30306 invoked from network); 15 Aug 2021 14:03:26 -0000 Received: from mother.openwall.net (195.42.179.200) by inbox.vuxu.org with ESMTPUTF8; 15 Aug 2021 14:03:26 -0000 Received: (qmail 12079 invoked by uid 550); 15 Aug 2021 14:03:25 -0000 Mailing-List: contact musl-help@lists.openwall.com; run by ezmlm Precedence: bulk List-Post: List-Help: List-Unsubscribe: List-Subscribe: List-ID: Reply-To: musl@lists.openwall.com Received: (qmail 12056 invoked from network); 15 Aug 2021 14:03:24 -0000 Date: Sun, 15 Aug 2021 10:03:12 -0400 From: Rich Felker To: Damian McGuckin Cc: musl@lists.openwall.com, Szabolcs Nagy Message-ID: <20210815140310.GG13220@brightrain.aerifal.cx> References: <0C6AAAD55DA44C6189B2FF4F5FB2C3E7@H270> <20210810213455.GB37904@port70.net> <20210814234612.GH37904@port70.net> <367A4018B58A4E308E2A95404362CBFB@H270> MIME-Version: 1.0 Content-Type: text/plain; charset=us-ascii Content-Disposition: inline In-Reply-To: User-Agent: Mutt/1.5.21 (2010-09-15) Subject: Re: [musl] [PATCH #2] Properly simplified nextafter() On Sun, Aug 15, 2021 at 06:24:31PM +1000, Damian McGuckin wrote: > > Hi Stefan, > > On Sun, 15 Aug 2021, Stefan Kanthak wrote: > > >__attribute__((noinline)) > >double nextafter(double x, double y) > >{ > >union {double f; unsigned long long i;} ux={x}, uy={y}; > >unsigned long long ax, ay; > >int e; > > > >if (isnan(x) || isnan(y)) > > return x + y; > >if (ux.i == uy.i) > > return y; > >#ifdef PATCH > >ax = ux.i << 1; > >ay = uy.i << 1; > >#else > >ax = ux.i & -1ULL/2; > >ay = uy.i & -1ULL/2; > >#endif > >if (ax == 0) { > > if (ay == 0) > > return y; > > ux.i = (uy.i & 1ULL<<63) | 1; > >#ifdef PATCH > >} else if (ax < ay == (long long) ux.i < 0) > >#else > >} else if (ax > ay || ((ux.i ^ uy.i) & 1ULL<<63)) > >#endif > > ux.i--; > >else > > ux.i++; > >e = ux.i >> 52 & 0x7ff; > >/* raise overflow if ux.f is infinite and x is finite */ > >if (e == 0x7ff) > > FORCE_EVAL(x + x); > >/* raise underflow if ux.f is subnormal or zero */ > >if (e == 0) > > FORCE_EVAL(x*x + ux.f*ux.f); > >return ux.f; > >} > > Maybe I am missing something and my brain is in weekend-mode ... > > I did a quick check and ran the above code for some test cases: > > nextafter(-9.7500000000e+01, 3.5000000000e+01) = -9.7500000000e+01 > yourpatch(-9.7500000000e+01, 3.5000000000e+01) = -9.7500000000e+01 > > The error is 2.8421709430e-14 > > nextafter(-9.7500000000e+01, -3.5000000000e+01) = -9.7500000000e+01 > yourpatch(-9.7500000000e+01, -3.5000000000e+01) = -9.7500000000e+01 > > The error is 2.8421709430e-14 I don't follow; are you claiming Stefan's patch introduces an error here? The outputs you printed show the exact same behavior before and after but it's possible you printed them wrong. The nextafter function is bit-exact; it does not have floating point error (inexactness) unless the implementation is buggy. > nextafter(-inf, inf) = -1.7976931349e+308 Correct > yourpatch(-inf, inf) = -nan Incorrect > > This is against standard GLIB. glibc I assume you mean? In any case yes this looks like a bug in the patch. Rich