From owner-freebsd-numerics@FreeBSD.ORG Mon Jun 1 23:33:41 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 0375089E for ; Mon, 1 Jun 2015 23:33:41 +0000 (UTC) (envelope-from stephen@missouri.edu) Received: from mst-rip5-missouri-out.um.umsystem.edu (mst-rip5-missouri-out.um.umsystem.edu [198.209.50.135]) (using TLSv1 with cipher RC4-SHA (128/128 bits)) (Client CN "um-tip1.um.umsystem.edu", Issuer "InCommon RSA Server CA" (verified OK)) by mx1.freebsd.org (Postfix) with ESMTPS id AC7C71DFE for ; Mon, 1 Jun 2015 23:33:40 +0000 (UTC) (envelope-from stephen@missouri.edu) X-IronPort-Anti-Spam-Filtered: true X-IronPort-Anti-Spam-Result: A2F+BgCK6mxV/8aeoM9cgxBUZIMYulhnAYEHVoYZgR09DwEBAQEBAQEDgQeCJSw7CAQJAQESAg0NWCMRVwEiAgUhAgQNAiEVEgQNCAEBiCmlA49fnz8BhF8ggSGSF4FFBZMQhDaISZVPI4N4gjWBAQEBAQ X-IPAS-Result: A2F+BgCK6mxV/8aeoM9cgxBUZIMYulhnAYEHVoYZgR09DwEBAQEBAQEDgQeCJSw7CAQJAQESAg0NWCMRVwEiAgUhAgQNAiEVEgQNCAEBiCmlA49fnz8BhF8ggSGSF4FFBZMQhDaISZVPI4N4gjWBAQEBAQ Received: from um-tcas2.um.umsystem.edu ([207.160.158.198]) by mst-rip5-exch-relay.um.umsystem.edu with ESMTP; 01 Jun 2015 18:32:31 -0500 Received: from UM-MBX-N02.um.umsystem.edu ([169.254.5.63]) by UM-TCAS2.um.umsystem.edu ([207.160.158.198]) with mapi id 14.03.0235.001; Mon, 1 Jun 2015 18:32:29 -0500 From: "Montgomery-Smith, Stephen" To: "freebsd-numerics@freebsd.org" Subject: Let's get moving Thread-Topic: Let's get moving Thread-Index: AQHQnMM5zb8ub/wH7U69fU25sa2Taw== Date: Mon, 1 Jun 2015 23:32:29 +0000 Message-ID: <556CEB8C.90605@missouri.edu> Accept-Language: en-US Content-Language: en-US X-MS-Has-Attach: X-MS-TNEF-Correlator: user-agent: Mozilla/5.0 (X11; Linux i686; rv:31.0) Gecko/20100101 Thunderbird/31.7.0 x-originating-ip: [207.160.158.194] Content-Type: text/plain; charset="utf-8" Content-ID: <9B76227003D1F84185B8E427E4FD40B0@missouri.edu> Content-Transfer-Encoding: base64 MIME-Version: 1.0 X-BeenThere: freebsd-numerics@freebsd.org X-Mailman-Version: 2.1.20 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: Mon, 01 Jun 2015 23:33:41 -0000 SGV5IHBlb3BsZSwNCg0KSSB3YXMgbG9va2luZyBhdCBodHRwczovL3dpa2kuZnJlZWJzZC5vcmcv TnVtZXJpY3MsIGFuZCB3ZSBhcmUgc3RhbGxpbmcuDQogRm9yIGV4YW1wbGUsIHdlIGhhdmUgd29y a2luZyBjb2RlIGZvciBjbG9nIGFuZCBjbG9nZiAtIHdoeSBkb24ndCB3ZQ0KY29tbWl0IGl0Pw0K DQpPbmUgb2YgdGhlIHJlYXNvbnMgSSBicmluZyB0aGlzIHVwIGlzIEkgd291bGQgbGlrZSB0byBj cmVhdGUgYSBwb3J0IGZvcg0KaHR0cDovL2xpYnJzYi5zb3VyY2Vmb3JnZS5uZXQvLCBhbmQgdGhl IG9ubHkgYmFycmllcnMgYXJlIHRoZSBsYWNrIG9mDQpjcG93ZiBhbmQgY3Bvdy4gIFJpZ2h0IG5v dyBJIGludGVuZCB0byBjb21taXQgaXQgd2l0aCBhIHJhdGhlcg0Kc2ltcGxpc3RpYyBwYXRjaDoN Cg0KK2NvbXBsZXggZmxvYXQgY2xvZ2YoY29tcGxleCBmbG9hdCBhKSB7DQorICAgICAgIHJldHVy biBsb2dmKGNhYnNmKGEpKSArIEkqY2FyZ2YoYSk7DQorfQ0KKw0KK2NvbXBsZXggZmxvYXQgY3Bv d2YoY29tcGxleCBmbG9hdCBhLCBjb21wbGV4IGZsb2F0IGIpIHsNCisgICAgICAgcmV0dXJuIGNl eHBmKGIqY2xvZ2YoYSkpOw0KK30NCisNCitjb21wbGV4IGRvdWJsZSBjbG9nKGNvbXBsZXggZG91 YmxlIGEpIHsNCisgICAgICAgcmV0dXJuIGxvZyhjYWJzKGEpKSArIEkqY2FyZyhhKTsNCit9DQor DQorY29tcGxleCBkb3VibGUgY3Bvdyhjb21wbGV4IGRvdWJsZSBhLCBjb21wbGV4IGRvdWJsZSBi KSB7DQorICAgICAgIHJldHVybiBjZXhwKGIqY2xvZyhhKSk7DQorfQ0KKw== From owner-freebsd-numerics@FreeBSD.ORG Tue Jun 2 00:40:37 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 522E1707 for ; Tue, 2 Jun 2015 00:40:37 +0000 (UTC) (envelope-from sgk@troutmask.apl.washington.edu) Received: from troutmask.apl.washington.edu (troutmask.apl.washington.edu [128.95.76.21]) (using TLSv1.2 with cipher ECDHE-RSA-AES256-GCM-SHA384 (256/256 bits)) (Client CN "troutmask", Issuer "troutmask" (not verified)) by mx1.freebsd.org (Postfix) with ESMTPS id 1730A1CD8 for ; Tue, 2 Jun 2015 00:40:37 +0000 (UTC) (envelope-from sgk@troutmask.apl.washington.edu) Received: from troutmask.apl.washington.edu (localhost [127.0.0.1]) by troutmask.apl.washington.edu (8.14.9/8.14.9) with ESMTP id t520QJjU045378 (version=TLSv1/SSLv3 cipher=DHE-RSA-AES256-GCM-SHA384 bits=256 verify=NO); Mon, 1 Jun 2015 17:26:20 -0700 (PDT) (envelope-from sgk@troutmask.apl.washington.edu) Received: (from sgk@localhost) by troutmask.apl.washington.edu (8.14.9/8.14.9/Submit) id t520QJNl045377; Mon, 1 Jun 2015 17:26:19 -0700 (PDT) (envelope-from sgk) Date: Mon, 1 Jun 2015 17:26:19 -0700 From: Steve Kargl To: "Montgomery-Smith, Stephen" Cc: "freebsd-numerics@freebsd.org" Subject: Re: Let's get moving Message-ID: <20150602002619.GA45253@troutmask.apl.washington.edu> References: <556CEB8C.90605@missouri.edu> MIME-Version: 1.0 Content-Type: text/plain; charset=us-ascii Content-Disposition: inline In-Reply-To: <556CEB8C.90605@missouri.edu> User-Agent: Mutt/1.5.23 (2014-03-12) X-BeenThere: freebsd-numerics@freebsd.org X-Mailman-Version: 2.1.20 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: Tue, 02 Jun 2015 00:40:37 -0000 On Mon, Jun 01, 2015 at 11:32:29PM +0000, Montgomery-Smith, Stephen wrote: > Hey people, > > I was looking at https://wiki.freebsd.org/Numerics, and we are stalling. The wiki is out-of-date. I've committed patches for erfcl, erfl, and lgammal. I'm working on powl and I'll do tgammal after I get powl done. I have versions of sincos[fl], which aren't quite ready, but I would like to have available for my work on the Bessel functions. > For example, we have working code for clog and clogf - why don't we > commit it? > > One of the reasons I bring this up is I would like to create a port for > http://librsb.sourceforge.net/, and the only barriers are the lack of > cpowf and cpow. Right now I intend to commit it with a rather > simplistic patch: >From src/math_private.h * The C99 standard intends x+I*y to be used for this, but x+I*y is * currently unusable in general since gcc introduces many overflow, * underflow, sign and efficiency bugs by rewriting I*y as * (0.0+I)*(y+0.0*I) and laboriously computing the full complex product. * In particular, I*Inf is corrupted to NaN+I*Inf, and I*-0 is corrupted * to -0.0+I*0.0. Last time I checked, clang was even worse with regards to I. -- Steve From owner-freebsd-numerics@FreeBSD.ORG Wed Jun 3 10:23:38 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 BD1AEAE1 for ; Wed, 3 Jun 2015 10:23:38 +0000 (UTC) (envelope-from brde@optusnet.com.au) Received: from mail106.syd.optusnet.com.au (mail106.syd.optusnet.com.au [211.29.132.42]) by mx1.freebsd.org (Postfix) with ESMTP id 8114F18EB for ; Wed, 3 Jun 2015 10:23:37 +0000 (UTC) (envelope-from brde@optusnet.com.au) Received: from c211-30-166-197.carlnfd1.nsw.optusnet.com.au (c211-30-166-197.carlnfd1.nsw.optusnet.com.au [211.30.166.197]) by mail106.syd.optusnet.com.au (Postfix) with ESMTPS id 38EDC3C4430; Wed, 3 Jun 2015 20:23:24 +1000 (AEST) Date: Wed, 3 Jun 2015 20:23:24 +1000 (EST) From: Bruce Evans X-X-Sender: bde@besplex.bde.org To: Steve Kargl cc: "Montgomery-Smith, Stephen" , "freebsd-numerics@freebsd.org" Subject: Re: Let's get moving In-Reply-To: <20150602002619.GA45253@troutmask.apl.washington.edu> Message-ID: <20150603201258.I1997@besplex.bde.org> References: <556CEB8C.90605@missouri.edu> <20150602002619.GA45253@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=XMDNMlVE c=1 sm=1 tr=0 a=KA6XNC2GZCFrdESI5ZmdjQ==:117 a=PO7r1zJSAAAA:8 a=kj9zAlcOel0A:10 a=JzwRw_2MAAAA:8 a=98Bt42YmqomoJD9UxmwA:9 a=CjuIK1q_8ugA:10 X-BeenThere: freebsd-numerics@freebsd.org X-Mailman-Version: 2.1.20 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: Wed, 03 Jun 2015 10:23:38 -0000 On Mon, 1 Jun 2015, Steve Kargl wrote: >> From src/math_private.h > > * The C99 standard intends x+I*y to be used for this, but x+I*y is > * currently unusable in general since gcc introduces many overflow, > * underflow, sign and efficiency bugs by rewriting I*y as > * (0.0+I)*(y+0.0*I) and laboriously computing the full complex product. > * In particular, I*Inf is corrupted to NaN+I*Inf, and I*-0 is corrupted > * to -0.0+I*0.0. > > Last time I checked, clang was even worse with regards to I. C11 now has CPLX() that can be used to handle this, and is used in a way that I don't like (to replace our old method). clang has some library support to do complex operations better. There are well known overflow problems for division, and even for addition if you want to let the real and imaginary parts overflow independently like we do in catrig*.c. The library functions might handle Infs and NaNs too. I haven't checked what they actually do. If you actually use them for every operation, then every operation becomes slow. libm wouldn't notice much since it mostly reduces to real and imaginary parts internally. Bruce From owner-freebsd-numerics@FreeBSD.ORG Wed Jun 3 10:29:06 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 A7E4AB6C for ; Wed, 3 Jun 2015 10:29:06 +0000 (UTC) (envelope-from brde@optusnet.com.au) Received: from mail108.syd.optusnet.com.au (mail108.syd.optusnet.com.au [211.29.132.59]) by mx1.freebsd.org (Postfix) with ESMTP id 399F01924 for ; Wed, 3 Jun 2015 10:29:05 +0000 (UTC) (envelope-from brde@optusnet.com.au) Received: from c211-30-166-197.carlnfd1.nsw.optusnet.com.au (c211-30-166-197.carlnfd1.nsw.optusnet.com.au [211.30.166.197]) by mail108.syd.optusnet.com.au (Postfix) with ESMTPS id DA8B31A1F66; Wed, 3 Jun 2015 20:12:08 +1000 (AEST) Date: Wed, 3 Jun 2015 20:12:05 +1000 (EST) From: Bruce Evans X-X-Sender: bde@besplex.bde.org To: "Montgomery-Smith, Stephen" cc: "freebsd-numerics@freebsd.org" Subject: Re: Let's get moving In-Reply-To: <556CEB8C.90605@missouri.edu> Message-ID: <20150603192222.Q1997@besplex.bde.org> References: <556CEB8C.90605@missouri.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=eZjABOwH c=1 sm=1 tr=0 a=KA6XNC2GZCFrdESI5ZmdjQ==:117 a=PO7r1zJSAAAA:8 a=kj9zAlcOel0A:10 a=JzwRw_2MAAAA:8 a=6I5d2MoRAAAA:8 a=FP58Ms26AAAA:8 a=kFOcfoAizQpEnSmb-jkA:9 a=CjuIK1q_8ugA:10 X-BeenThere: freebsd-numerics@freebsd.org X-Mailman-Version: 2.1.20 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: Wed, 03 Jun 2015 10:29:06 -0000 On Mon, 1 Jun 2015, Montgomery-Smith, Stephen wrote: > Hey people, > > I was looking at https://wiki.freebsd.org/Numerics, and we are stalling. > For example, we have working code for clog and clogf - why don't we > commit it? I thought that you might commit these before your complex inverse hyperbolic functions :-). catrig.c and catrigf.c have been committed. These pass my accuracy tests, but I'm not happy with their integration. Bugs start in their file names. All other files in msun/src have names with prefixes like e_ or s_. Comments are stripped in the committed version of catrigf.c, and your template for maintaining the comments was not committed. catrigl.c was not committed. I only have incomplete minor patches for this. Most were in patches that I sent you in 2012 to merge: X --- stephen-catrigl.c 2012-09-22 21:14:24.000000000 +0000 X +++ catrigl.c 2013-06-08 10:47:56.000000000 +0000 X @@ -60,6 +69,6 @@ X #if LDBL_MANT_DIG == 64 X static const union IEEEl2bits X -um_e = LD80C(0xadf85458a2bb4a9b, 1, 0, 2.71828182845904523536e0L), X -um_ln2 = LD80C(0xb17217f7d1cf79ac, -1, 0, 6.93147180559945309417e-1L); X +um_e = LD80C(0xadf85458a2bb4a9b, 1, 2.71828182845904523536e+0L), X +um_ln2 = LD80C(0xb17217f7d1cf79ac, -1, 6.93147180559945309417e-1L); X #define m_e um_e.e X #define m_ln2 um_ln2.e Catch up with API changes in LD80C(). These have been in -current since 2012 so are standard now. X @@ -348,5 +357,5 @@ X X if (y == 0 && ax <= 1) X - return (cpackl(atanhl(x), y)); /* XXX need atanhl() */ X + return (cpackl(atanhl(x), y)); X X if (x == 0) In 2012, this used our hacked up versions of atanhl(). A production version was committed in 2013. Someone churned the API in math_private.h so that cpackl() no longer exists. I don't used the churned parts, so the above still compiles for me. X @@ -369,12 +378,7 @@ X } X X - if (ax == 1 && ay < LDBL_EPSILON) { X -#if 0 X - if (ay > 2*LDBL_MIN) X - rx = - logl(ay/2) / 2; X - else X -#endif X - rx = - (logl(ay) - m_ln2) / 2; X - } else X + if (ax == 1 && ay < LDBL_EPSILON) X + rx = (logl(ay) - m_ln2) / -2; X + else X rx = log1pl(4*ax / sum_squares(ax-1, ay)) / 4; X Minor cleanups. IIRC, this was not in the template, and the float, double and long versions had difference here from being manually edited. The committed float and double versions have this and some other cleanups, especially in the float version where generation from the template made a mess. The main things missing here are: - no support for i387 mode switching (ENTERI()). Not very easy to add cleanly since there are many entry and exit points. - not up to date with cpackl(). My version of clog*() is a little more finished. It has the i387 mode switching support for clogl(), but is not up to date with cpack*(). It needs splitting into separate files for each precision. It handles comment duplication not very well by just duplicating comments and manually maintaining consistency. catrig*() can also use clog*() when it exists, but have fairly good private versions for the special cases needed. > One of the reasons I bring this up is I would like to create a port for > http://librsb.sourceforge.net/, and the only barriers are the lack of > cpowf and cpow. Right now I intend to commit it with a rather > simplistic patch: > > +complex float clogf(complex float a) { > + return logf(cabsf(a)) + I*cargf(a); > +} > + > +complex float cpowf(complex float a, complex float b) { > + return cexpf(b*clogf(a)); > +} > + > +complex double clog(complex double a) { > + return log(cabs(a)) + I*carg(a); > +} > + > +complex double cpow(complex double a, complex double b) { > + return cexp(b*clog(a)); > +} I think a lot of ports already have similar hacks. libm has similar hacks for powl() and tgammal(), to support C++. It doesn't bother with similar hacks for complex functions. C99 and at least old POSIX don't have complex Bessel or gamma functions, so we're closer to having all standard real functions that complex ones. Only powl() and cpowl() are standard, unsupported, not written, and difficult. Bruce