From owner-freebsd-numerics@FreeBSD.ORG Fri Jan 30 05:07:04 2015 Return-Path: Delivered-To: freebsd-numerics@FreeBSD.org Received: from mx1.freebsd.org (mx1.freebsd.org [IPv6:2001:1900:2254:206a::19:1]) (using TLSv1.2 with cipher AECDH-AES256-SHA (256/256 bits)) (No client certificate requested) by hub.freebsd.org (Postfix) with ESMTPS id 28C79EA for ; Fri, 30 Jan 2015 05:07:04 +0000 (UTC) Received: from mail.apache.org (hermes.apache.org [140.211.11.3]) by mx1.freebsd.org (Postfix) with SMTP id 07444F4D for ; Fri, 30 Jan 2015 05:07:03 +0000 (UTC) Received: (qmail 74098 invoked by uid 99); 30 Jan 2015 05:07:02 -0000 Received: from mail-relay.apache.org (HELO mail-relay.apache.org) (140.211.11.15) by apache.org (qpsmtpd/0.29) with ESMTP; Fri, 30 Jan 2015 05:07:02 +0000 Received: from [192.168.0.103] (unknown [190.157.136.22]) by mail-relay.apache.org (ASF Mail Server at mail-relay.apache.org) with ESMTPSA id C21961A0386 for ; Fri, 30 Jan 2015 05:07:01 +0000 (UTC) From: Pedro Giffuni Subject: Something wrong in scalblnl() Message-Id: Date: Fri, 30 Jan 2015 00:06:58 -0500 To: freebsd-numerics@FreeBSD.org Mime-Version: 1.0 (Mac OS X Mail 8.2 \(2070.6\)) X-Mailer: Apple Mail (2.2070.6) Content-Type: text/plain; charset=utf-8 Content-Transfer-Encoding: quoted-printable X-Content-Filtered-By: Mailman/MimeDel 2.1.18-1 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: Fri, 30 Jan 2015 05:07:04 -0000 Dear FreeBSD math gurus; Looking at the clang checker report: = http://scan.freebsd.org/scan-build/lib/2015-01-05-amd64/report-4c1ee4.html= #EndPath = in the scalblnl() implementation we are assigning a value to =E2=80=9Cin=E2= =80=9D and but we are ignoring it for the calculation. looking at the = family of similar functions nearby, it looks clear that we have a typo = in the return value. Pedro.= From owner-freebsd-numerics@FreeBSD.ORG Fri Jan 30 07:40:00 2015 Return-Path: Delivered-To: freebsd-numerics@freebsd.org Received: from mx1.freebsd.org (mx1.freebsd.org [8.8.178.115]) (using TLSv1.2 with cipher AECDH-AES256-SHA (256/256 bits)) (No client certificate requested) by hub.freebsd.org (Postfix) with ESMTPS id 8BD0265B; Fri, 30 Jan 2015 07:40:00 +0000 (UTC) Received: from mail108.syd.optusnet.com.au (mail108.syd.optusnet.com.au [211.29.132.59]) by mx1.freebsd.org (Postfix) with ESMTP id 4DEDEF21; Fri, 30 Jan 2015 07:39:59 +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 mail108.syd.optusnet.com.au (Postfix) with ESMTPS id 79B6B1A2D94; Fri, 30 Jan 2015 18:39:48 +1100 (AEDT) Date: Fri, 30 Jan 2015 18:39:48 +1100 (EST) From: Bruce Evans X-X-Sender: bde@besplex.bde.org To: Pedro Giffuni Subject: Re: Something wrong in scalblnl() In-Reply-To: Message-ID: <20150130182211.D1663@besplex.bde.org> References: MIME-Version: 1.0 X-Optus-CM-Score: 0 X-Optus-CM-Analysis: v=2.1 cv=E/j0+8tl c=1 sm=1 tr=0 a=7NqvjVvQucbO2RlWB8PEog==:117 a=PO7r1zJSAAAA:8 a=JzwRw_2MAAAA:8 a=nlC_4_pT8q9DhB4Ho9EA:9 a=cz2ZRIgtxKwA:10 a=wJWlkF7cXJYA:10 a=c3-DdYJoA5YA:10 a=6I5d2MoRAAAA:8 a=p-n1nNFXNEepaHI2oUoA:9 a=45ClL6m2LaAA:10 a=QhiZzMmKTVwA:10 Content-Type: TEXT/PLAIN; charset=X-UNKNOWN; format=flowed Content-Transfer-Encoding: QUOTED-PRINTABLE X-Content-Filtered-By: Mailman/MimeDel 2.1.18-1 Cc: freebsd-numerics@freebsd.org 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: Fri, 30 Jan 2015 07:40:00 -0000 On Fri, 30 Jan 2015, Pedro Giffuni wrote: > Looking at the clang checker report: > > http://scan.freebsd.org/scan-build/lib/2015-01-05-amd64/report-4c1ee4.htm= l#EndPath > > in the scalblnl() implementation we are assigning a value to =E2=80=9Cin= =E2=80=9D and but we are ignoring it for the calculation. looking at the fa= mily of similar functions nearby, it looks clear that we have a typo in the= return value. The second-last line should use 'in', as in the other functions. This shows that exceptional cases in scalbnl() were never tested. Things to fix nearby: - gnu-style parentheses - the conversion can be done a bit more cleanly using =09in =3D (n > MAX) ? MAX : (n < MIN) ? MIN : n; where MAX and MIN can be INT_MAX and INT_MIN or just +- of a large value related to the maximum exponent. Bruce From owner-freebsd-numerics@FreeBSD.ORG Fri Jan 30 22:39:36 2015 Return-Path: Delivered-To: freebsd-numerics@freebsd.org Received: from mx1.freebsd.org (mx1.freebsd.org [8.8.178.115]) (using TLSv1.2 with cipher AECDH-AES256-SHA (256/256 bits)) (No client certificate requested) by hub.freebsd.org (Postfix) with ESMTPS id 1A318A51; Fri, 30 Jan 2015 22:39:36 +0000 (UTC) Received: from vps.rulingia.com (vps.rulingia.com [103.243.244.15]) (using TLSv1.2 with cipher ECDHE-RSA-AES256-GCM-SHA384 (256/256 bits)) (Client CN "vps.rulingia.com", Issuer "CAcert Class 3 Root" (not verified)) by mx1.freebsd.org (Postfix) with ESMTPS id A8372100; Fri, 30 Jan 2015 22:39:35 +0000 (UTC) Received: from server.rulingia.com (c220-239-242-83.belrs5.nsw.optusnet.com.au [220.239.242.83]) by vps.rulingia.com (8.14.9/8.14.9) with ESMTP id t0UMdHAK079701 (version=TLSv1/SSLv3 cipher=DHE-RSA-AES256-GCM-SHA384 bits=256 verify=OK); Sat, 31 Jan 2015 09:39:23 +1100 (AEDT) (envelope-from peter@rulingia.com) X-Bogosity: Ham, spamicity=0.000000 Received: from server.rulingia.com (localhost.rulingia.com [127.0.0.1]) by server.rulingia.com (8.14.9/8.14.9) with ESMTP id t0UMdAOC010935 (version=TLSv1/SSLv3 cipher=DHE-RSA-AES256-GCM-SHA384 bits=256 verify=NO); Sat, 31 Jan 2015 09:39:10 +1100 (AEDT) (envelope-from peter@server.rulingia.com) Received: (from peter@localhost) by server.rulingia.com (8.14.9/8.14.9/Submit) id t0UMd9Ed010934; Sat, 31 Jan 2015 09:39:09 +1100 (AEDT) (envelope-from peter) Date: Sat, 31 Jan 2015 09:39:09 +1100 From: Peter Jeremy To: Steve Kargl Subject: Re: clang is almost useless for complex arithmetic Message-ID: <20150130223909.GA10817@server.rulingia.com> References: <20140326002205.GA9940@troutmask.apl.washington.edu> MIME-Version: 1.0 Content-Type: multipart/signed; micalg=pgp-sha512; protocol="application/pgp-signature"; boundary="zhXaljGHf11kAtnf" Content-Disposition: inline In-Reply-To: <20140326002205.GA9940@troutmask.apl.washington.edu> X-PGP-Key: http://www.rulingia.com/keys/peter.pgp User-Agent: Mutt/1.5.23 (2014-03-12) Cc: freebsd-numerics@freebsd.org, freebsd-toolchain@freebsd.org 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: Fri, 30 Jan 2015 22:39:36 -0000 --zhXaljGHf11kAtnf Content-Type: text/plain; charset=us-ascii Content-Disposition: inline Content-Transfer-Encoding: quoted-printable To resurrect an old thread... On 2014-Mar-25 17:22:05 -0700, Steve Kargl wrote: >It appears that clang developers have chosen the naive >complex division algorithm, and it does not matter whether >one turns CX_LIMITED_RANGE on or off. This means that=20 >if one uses clang with complex types, one must be careful >with the range of values allowed in complex division. In >other words, implementation of complex libm routines cannot >use complex data types and must fallback to a decomposition >into real and imaginary components.=20 Whilst the fixes don't seem to have made it into FreeBSD yet, it seems that this has been improved in recent Clang/LLVM - see http://reviews.llvm.org/D5698, committed as rL219557. I re-implemented their division algorithm in C and checked it against the algorithms listed in http://arxiv.org/pdf/1210.4539v2.pdf. Whilst it fails on same of the difficult values listed in that paper, when faced with a range of random arguments, it seems to perform better[*] than the "robust" algorithm. I didn't look closely at the multiplication algorithm but, based on a quick look, it seems that it's still susceptable to spurious overflows. [*] or at least closer to the result given by the naive algorithm using x87 long doubles. --=20 Peter Jeremy --zhXaljGHf11kAtnf Content-Type: application/pgp-signature -----BEGIN PGP SIGNATURE----- Version: GnuPG v2 iQJ8BAEBCgBmBQJUzAgNXxSAAAAAAC4AKGlzc3Vlci1mcHJAbm90YXRpb25zLm9w ZW5wZ3AuZmlmdGhob3JzZW1hbi5uZXRFRUIyOTg2QzMwNjcxRTc0RTY1QzIyN0Ux NkE1OTdBMEU0QTIwQjM0AAoJEBall6Dkogs0NHwP/RRvnj9dVpunArX8WNe/eTty JWtsL0OiQCmRxvUa+1nmAggUJxxXxjEhwx2CwrBRZEBppKOdxrCrcUxd6lwa56Xz d/phocOUckfbDT4JaiksQHhucba0fXW06MuAiqnc6pYEOvWj715Lz3VkICJG/ZMZ sIq6R0g+/4ZFXMHybVBH65hiiIXzFiWIIyx2XyOb14zrgCxPHPGpQDPiddiuQTPD cB1WLcNUpKDwAR0f8WFsvZLLxJ8Zj23IyjMNzI2Tlb1IhoCTRUmMYdXZyg3ai1XY fqpts+pR2JMTDddR4Sp1W1uwPKYASFx7wxZOkWyd52YMx7zvvNHuPhc3HAvVxJxw +cePrV6ArlWU+/pEXO5KlsYtB2XLbqLLSzsKEVOWC8eOUWmjs0aTPXqHQ3SLzoGz J5y4K6rNc1z+kuHBOoJItYXAL8GYLUDC13ZnDX7DYTfTZHE5JN2sKs5FZPO8GujS hGzCm3D6xvk6+zbQcrTzS5C+wQtYc5kKfM2EqcHfcY0R0kMWFYO5qe4mia2oP+xk 4uSahOYPmnuLM3wrOrOO+y+GBLIdkJPEZc8w3UZtxhULC3h3ohHLKdfQQGBvyeG9 00o0jNH4hJ0SFT9u6gX2EfHv0nEwlwZzzaPsvTpQy9TyoVsGY4ipduheBqzaZ+kW 1IaujlvmTTAfcqP8Aoms =VeVJ -----END PGP SIGNATURE----- --zhXaljGHf11kAtnf-- From owner-freebsd-numerics@FreeBSD.ORG Sat Jan 31 08:32:11 2015 Return-Path: Delivered-To: freebsd-numerics@freebsd.org Received: from mx1.freebsd.org (mx1.freebsd.org [IPv6:2001:1900:2254:206a::19:1]) (using TLSv1.2 with cipher AECDH-AES256-SHA (256/256 bits)) (No client certificate requested) by hub.freebsd.org (Postfix) with ESMTPS id A0075E7F; Sat, 31 Jan 2015 08:32:11 +0000 (UTC) Received: from mail104.syd.optusnet.com.au (mail104.syd.optusnet.com.au [211.29.132.246]) by mx1.freebsd.org (Postfix) with ESMTP id 61A4FE84; Sat, 31 Jan 2015 08:32:11 +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 mail104.syd.optusnet.com.au (Postfix) with ESMTPS id CF8CC426E99; Sat, 31 Jan 2015 19:32:01 +1100 (AEDT) Date: Sat, 31 Jan 2015 19:31:58 +1100 (EST) From: Bruce Evans X-X-Sender: bde@besplex.bde.org To: Peter Jeremy Subject: Re: clang is almost useless for complex arithmetic In-Reply-To: <20150130223909.GA10817@server.rulingia.com> Message-ID: <20150131184056.M1858@besplex.bde.org> References: <20140326002205.GA9940@troutmask.apl.washington.edu> <20150130223909.GA10817@server.rulingia.com> 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=baJSDo/B c=1 sm=1 tr=0 a=7NqvjVvQucbO2RlWB8PEog==:117 a=PO7r1zJSAAAA:8 a=kj9zAlcOel0A:10 a=JzwRw_2MAAAA:8 a=Twlkf-z8AAAA:8 a=EWPDJS0nAAAA:8 a=eOTlJJCRNr-8Tytd1pgA:9 a=CjuIK1q_8ugA:10 Cc: freebsd-numerics@freebsd.org, Steve Kargl , freebsd-toolchain@freebsd.org 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, 31 Jan 2015 08:32:11 -0000 On Sat, 31 Jan 2015, Peter Jeremy wrote: > To resurrect an old thread... > > On 2014-Mar-25 17:22:05 -0700, Steve Kargl wrote: >> It appears that clang developers have chosen the naive >> complex division algorithm, and it does not matter whether >> one turns CX_LIMITED_RANGE on or off. This means that >> if one uses clang with complex types, one must be careful >> with the range of values allowed in complex division. In >> other words, implementation of complex libm routines cannot >> use complex data types and must fallback to a decomposition >> into real and imaginary components. > > Whilst the fixes don't seem to have made it into FreeBSD yet, it seems > that this has been improved in recent Clang/LLVM - see > http://reviews.llvm.org/D5698, committed as rL219557. The libcalls (as in gcc) allow replacing the method easily. > I re-implemented their division algorithm in C and checked it against > the algorithms listed in http://arxiv.org/pdf/1210.4539v2.pdf. Whilst > it fails on same of the difficult values listed in that paper, when > faced with a range of random arguments, it seems to perform better[*] > than the "robust" algorithm. > > I didn't look closely at the multiplication algorithm but, based on a > quick look, it seems that it's still susceptable to spurious overflows. > > [*] or at least closer to the result given by the naive algorithm using > x87 long doubles. Using long doubles would fix many cases without even trying. But clang breaks their automatic use even on i386 with -march optimizations, by using SSE instructions if the arch supports them. With gcc, it takes more magic compiler options (-mfpmath=sse+387) to get this bug and this is still documented as experimental. The libcalls could use long doubles even more intentionally, but then in FreeBSD they would have to fight with the default rounding precision being double (so long doubles are much like doubles and using them only helps for float complex). I checked some newer gcc (4.8) support for extra precision. There is a -fextra-precision=style option, except it is spelled -fexcess-precision. -fexcess-precision=fast gives the historical broken behaviour. -fexcess-precision=standard is claimed to give the C99 behaviour. It may give it, but it also gives the C11 bug of destroying extra precision on function return (including for inline functions). -std=c99 implies -fexcess-precision=standard. This is good for purity, but mostly bad. It is too surprising that -std=c99, gives a large pessimization. Fortunately, -std=c99 is often unusuable for other reasons. FreeBSD uses -std=gnu99 and that gives the historical behaviour. -fexcess-precision=standard is several years old now, but not supported by clang. Destruction of extra precision on return is not done by clang even with -std=c11. This goes with not supporting C99 yet. Bruce