[GHC] #14927: Hyperbolic area sine is unstable for (even moderately) big negative arguments.

GHC ghc-devs at haskell.org
Thu Mar 15 17:43:50 UTC 2018


#14927: Hyperbolic area sine is unstable for (even moderately) big negative
arguments.
-------------------------------------+-------------------------------------
           Reporter:                 |             Owner:  (none)
  leftaroundabout                    |
               Type:  bug            |            Status:  new
           Priority:  normal         |         Milestone:
          Component:                 |           Version:  8.2.1
  libraries/base                     |
           Keywords:                 |  Operating System:  Unknown/Multiple
       Architecture:                 |   Type of failure:  Incorrect result
  Unknown/Multiple                   |  at runtime
          Test Case:                 |        Blocked By:
           Blocking:                 |   Related Tickets:
Differential Rev(s):                 |         Wiki Page:
-------------------------------------+-------------------------------------
 `asinh` is supposed to be the inverse of `sinh`, and this works pretty
 reliable for positive arguments. However, for negative arguments, the
 [http://hackage.haskell.org/package/base-4.10.1.0/docs/src/GHC.Float.html#line-473
 currently used formula]
 {{{
 asinh x = log (x + sqrt (1.0+x*x))
 }}}
 gets unstable much earlier than necessary, namely when the summands in the
 logarithm cancel almost to zero, dominated by the numerical error of the
 square root.

 This is particularly troubling because mathematically **a)** `asinh` is a
 very “inert” function (i.e. you can carelessly put in huge numbers and –
 as long as they're not outright `Infinity` – always get a somewhat sane
 result), pseudo-sigmoidal as it were **b)** it is an ''odd function'',
 fulfilling `asinh (-x) = -asinh x`.

 Both is reflected in other implementations, e.g. Python, but not in GHC
 Haskell:
 {{{
 GHCi, version 8.2.1         Python 3.5.2 (default, Nov 23 2017, 16:37:01)

                             In [1]: from math import *

 Prelude> asinh 1e6          In [2]: asinh(1e6)
 14.50865773852447           Out[2]: 14.50865773852447

 Prelude> asinh (-1e6)       In [3]: asinh(-1e6)
 -14.50865012405984          Out[3]: -14.50865773852447

 Prelude> asinh 1e9          In [4]: asinh(1e9)
 21.416413017506358          Out[4]: 21.416413017506354

 Prelude> asinh (-1e9)       In [5]: asinh(-1e9)
 -Infinity                   Out[5]: -21.416413017506354

 Prelude> asinh 1e76         In [6]: asinh(1e76)
 175.6896142481074           Out[6]: 175.68961424810743

 Prelude> asinh (-1e76)      In [7]: asinh(-1e76)
 -Infinity                   Out[7]: -175.68961424810743
 }}}

 Demo of non-inverse property:
 {{{
 Prelude> [(x, asinh $ sinh x) | x <- [-25..25]]
 [(-25.0,-Infinity)
 ,(-24.0,-Infinity)
 ,(-23.0,-Infinity)
 ,(-22.0,-Infinity)
 ,(-21.0,-Infinity)
 ,(-20.0,-Infinity)
 ,(-19.0,-18.021826694558577)
 ,(-18.0,-18.021826694558577)
 ,(-17.0,-17.0102257828801)
 ,(-16.0,-15.998624871201619)
 ,(-15.0,-14.999878578873695)
 ,(-14.0,-13.999968823323222)
 ,(-13.0,-12.999991335176079)
 ,(-12.0,-12.000000137072186)
 ,(-11.0,-10.999999903206444)
 ,(-10.0,-10.000000013503529)
 ,(-9.0,-9.000000000551713)
 ,(-8.0,-8.00000000017109)
 ,(-7.0,-7.000000000036329)
 ,(-6.0,-5.999999999998066)
 ,(-5.0,-5.000000000000641)
 ,(-4.0,-4.000000000000046)
 ,(-3.0,-2.999999999999989)
 ,(-2.0,-1.9999999999999991)
 ,(-1.0,-1.0)
 ,(0.0,0.0)
 ,(1.0,1.0)
 ,(2.0,2.0)
 ,(3.0,3.0)
 ,(4.0,4.0)
 ,(5.0,5.0)
 ,(6.0,6.0)
 ,(7.0,7.0)
 ,(8.0,8.0)
 ,(9.0,9.0)
 ,(10.0,10.0)
 ,(11.0,11.0)
 ,(12.0,12.0)
 ,(13.0,13.0)
 ,(14.0,14.0)
 ,(15.0,15.0)
 ,(16.0,16.0)
 ,(17.0,17.0)
 ,(18.0,18.0)
 ,(19.0,19.0)
 ,(20.0,20.0)
 ,(21.0,21.0)
 ,(22.0,22.0)
 ,(23.0,23.0)
 ,(24.0,24.0)
 ,(25.0,25.0)]
 }}}
 Those results are less than satisfying, even for inputs that aren't
 astronomically big at all.

 A simple fix would be to “copy” the sane positive-number behaviour to the
 negative side:
 {{{
     asinh x
       | x < 0      = -asinh (-x)
       | otherwise  = log (x + sqrt (1.0+x*x))
 }}}

-- 
Ticket URL: <http://ghc.haskell.org/trac/ghc/ticket/14927>
GHC <http://www.haskell.org/ghc/>
The Glasgow Haskell Compiler


More information about the ghc-tickets mailing list