accuracy of asinh and atanh

Previous Topic Next Topic
 
classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

accuracy of asinh and atanh

Matt Peddie
Hi devs,

I tried to use  asinh :: Double -> Double  and discovered that it's
inaccurate compared to my system library (GNU libm), even returning
-Infinity in place of finite values in the neighborhood of -22 for
large negative arguments.  `atanh` is also inaccurate compared to the
system library.  I wrote up a more detailed description of the problem
including plots in the README file at
https://github.com/peddie/ghc-inverse-hyperbolic -- this repository is
package that can help you examine the error for yourself or generate
the plots, and it also contains accurate pure-Haskell translations of
the system library's implementation for these functions.  What's the
next step to fixing this in GHC?

Cheers

Matt Peddie
_______________________________________________
ghc-devs mailing list
[hidden email]
http://mail.haskell.org/cgi-bin/mailman/listinfo/ghc-devs
Reply | Threaded
Open this post in threaded view
|

Re: accuracy of asinh and atanh

Edward Kmett-2
Note: From skimming your readme it is worth noting that  log1p _is_ in base now (alongside expm1, log1pexp, and log1mexp). We added them all a couple of years back as a result of the very thread linked in your README.

You need to `import Numeric` to see them, though.

Switching to more accurate functions for doubles and floats for asinh, atanh, etc. to exploit this sort of functionality at least seems to make a lot of sense.

That can be done locally without any user API impact as the current definitions aren't supplied as defaults, merely as pointwise implementations instance by instance. Things will just become more accurate.

In that same spirit, we can probably crib a better version for complex numbers from somewhere as well, as it follows the same general simplistic formula right now, even if it can't be plugged directly into the equations you've given. For that matter, the log1p definition we're using for complex numbers was the best I could come up with, but there may well be a more accurate version you can find down in the mines of libm or another math library written by real analysts.
    log1p x@(a :+ b)
      | abs a < 0.5 && abs b < 0.5
      , u <- 2*a + a*a + b*b = log1p (u/(1 + sqrt(u+1))) :+ atan2 (1 + a) b
      | otherwise = log (1 + x)

So, here's a +1 from the libraries committee side towards improving the situation.

From there, it's a small matter of implementation. 

Here's where I'd usually get Ben involved. Hi Ben!

-Edward

On Sat, Jun 2, 2018 at 1:23 AM, Matt Peddie <[hidden email]> wrote:
Hi devs,

I tried to use  asinh :: Double -> Double  and discovered that it's
inaccurate compared to my system library (GNU libm), even returning
-Infinity in place of finite values in the neighborhood of -22 for
large negative arguments.  `atanh` is also inaccurate compared to the
system library.  I wrote up a more detailed description of the problem
including plots in the README file at
https://github.com/peddie/ghc-inverse-hyperbolic -- this repository is
package that can help you examine the error for yourself or generate
the plots, and it also contains accurate pure-Haskell translations of
the system library's implementation for these functions.  What's the
next step to fixing this in GHC?

Cheers

Matt Peddie
_______________________________________________
ghc-devs mailing list
[hidden email]
http://mail.haskell.org/cgi-bin/mailman/listinfo/ghc-devs


_______________________________________________
ghc-devs mailing list
[hidden email]
http://mail.haskell.org/cgi-bin/mailman/listinfo/ghc-devs
Reply | Threaded
Open this post in threaded view
|

Re: accuracy of asinh and atanh

Ben Gamari-3
Sending again, this time including ghc-devs.

Edward Kmett <[hidden email]> writes:

> Note: From skimming your readme it is worth noting that  log1p _is_ in base
> now (alongside expm1, log1pexp, and log1mexp). We added them all a couple
> of years back as a result of the very thread linked in your README.
>
> You need to `import Numeric` to see them, though.
>
> Switching to more accurate functions for doubles and floats for asinh,
> atanh, etc. to exploit this sort of functionality at least seems to make a
> lot of sense.
>
> That can be done locally without any user API impact as the current
> definitions aren't supplied as defaults, merely as pointwise
> implementations instance by instance. Things will just become more accurate.
>
> In that same spirit, we can probably crib a better version for complex
> numbers from somewhere as well, as it follows the same general simplistic
> formula right now, even if it can't be plugged directly into the equations
> you've given. For that matter, the log1p definition we're using for complex
> numbers was the best I could come up with, but there may well be a more
> accurate version you can find down in the mines of libm or another math
> library written by real analysts.
>
[snip]
>
> So, here's a +1 from the libraries committee side towards improving the
> situation.
>
> From there, it's a small matter of implementation.
>
> Here's where I'd usually get Ben involved. Hi Ben!
>
Hi Edward and Matt!

Indeed the sinh sinh situation should already be improved in 8.6. See
#14927 and GHC commit 3ea33411d7cbf32c20940cc72ca07df6830eeed7. I
suspect this should fix Matt's issue.

Regarding log1p, I'd happily accept patches improving on the status quo.

Cheers,

- Ben

_______________________________________________
ghc-devs mailing list
[hidden email]
http://mail.haskell.org/cgi-bin/mailman/listinfo/ghc-devs

signature.asc (668 bytes) Download Attachment