Use of C99 extra long double math functions after r236148
Steve Kargl
sgk at troutmask.apl.washington.edu
Sun Aug 12 22:57:45 UTC 2012
On Wed, Jul 18, 2012 at 08:53:28AM +1000, Peter Jeremy wrote:
> On 2012-Jul-16 21:38:48 -0700, Steve Kargl <sgk at troutmask.apl.washington.edu> wrote:
> >On Tue, Jul 17, 2012 at 02:21:25PM +1000, Peter Jeremy wrote:
> >> Again, handling the special cases listed in G.6 is all just
> >> boilerplate code that we can take as assumed for pseudocode. IMO, it
> >
> >As someone who spent 10+ years getting sqrtl(), cbrtl(), ccosh(), sinl(),
> >cosl(), tanl(), etc into FreeBSD, I respectfully disagree with your
> >take that it is just boilerplate. Getting this stuff right is much harder
> >than I think some people understand. Oh well, I'll back to lurking and
> >working on things I need.
>
> OK. I'll admit that I haven't tried this before but attached is my
> first try at catan[h](). It should cover all the special casing
> according to WG14/N1256 and compiles if you add the relevant
> declarations to complex.h. I'd appreciate feedback.
>
> Notes on it:
> - The actual code at the end of catanh() is just one possible
> algorithm. It's not intended as final code and the final code is
> likely to need additional special case handling to minimise
> precision loss and prevent unwanted exceptions.
> - cpack(-cimag(r), -creal(r)) gives better code than
> -cpack(cimag(r), creal(r)) on i386 and identical code on amd64.
> - The fpclassify() macros are bitmasks on FreeBSD. Assuming this
> would allow (ci == FP_ZERO || ci == FP_NAN) to be simplified but the
> standard only requires that they have distinct values and (eg)
> Solaris implements them as a series so the current code is more
> portable.
>
I won't have time to go over the code in detail until
this weekend, but a quick peek showed some issues. The
first is style. Although fdlibm has a rather interest
coding style, new code should use KNF.
/*
* Calculate complex arc tangent using the identity:
* catan(z) = -i catanh(iz)
*/
double complex
catan(double complex z)
{
complex double r;
r = catanh(cpack(cimag(z), creal(z)));
I think you're missing a sign. Let z = x + i*y.
Then, i*z = i*x+i*i*y = -y + i*x, yielding
r = catanh(cpack(-cimag(z), creal(z)));
return (cpack(-cimag(r), -creal(r)));
Again, it seems a sign error has occurred.
Let catanh(i*z) = u + i*v. Then, you have
-i*catanh(i*z) = -i*u-i*i*v = v-i*u, yielding
return (cpack(cimag(r), -creal(r)));
}
Of coure, I could be wrong.
--
Steve
More information about the freebsd-numerics
mailing list