While adding complex numbers to WebRacket I have found a potential problem
in the math/base implementation of atanh. When called with complex
arguments, the accuracy isn't as good as Chez Scheme.
For the complex case the implementation of atanh in math/base uses:
(* 0.5 (- (log (+ 1.0 x)) (log (- 1.0 x))))
This is mathematically correct, but not numerically stable.
The problem is how to improve the situation?
Option 1) improve the implementation in math/base.
Option 2) figure out how to call the Chez Scheme primitive atanh.
The section below was produced by Codex when I asked it examine the
accuracy of atanh on complex numbers.
Note: There might be similar problems with asinh and and acosh.
/Jens Axel
Codex:
For math/base, complex accuracy is not uniformly good. It uses the direct formula at [base-functions.rkt (line 104)](/Users/soegaard/Dropbox/GitHub/racket/racket/share/pkgs/math-lib/math/private/base/base-functions.rkt:104):
(* 0.5 (- (log (+ 1.0 x)) (log (- 1.0 x))))
That is mathematically right, but numerically fragile.
Measured against a 4096-bit bigfloat reference:
Region | math/base max-component error
-- | --
ordinary random [-10,10]^2 | median 5 ulp, p95 27 ulp, worst 23515 ulp
tiny complex inputs | often loses the real part entirely
near 1+iε | very good, usually 0-1 ulp
huge complex inputs | often loses the small real part entirely
Examples:
z = 1e-20+1e-20i
ref 1e-20+1e-20i
math/base 0.0+1e-20i
z = 1e100+0.0i
ref 1e-100+1.5707963267948966i
math/base 0.0+1.5707963267948966i
The Chez source-level cflatanh at 5_3.ss (line 761) is much more accuracy-oriented: it uses Kahan-style formulas and, when log1p is available through fllog1+ (line 412), it preserves those tiny real parts. In my model of that path, ordinary/tiny cases were typically 0-1 ulp.
But Chez has its own weak spot near z = 1+iε: it adds rho ≈ 2.98e-154 to |imag|, so for ε << rho it clamps the behavior:
z = 1.0+1e-300i
ref 345.73433753938684+0.7853981633974483i
math/base 345.7343375393868 +0.7853981633974483i
Chez model 177.09910463306602+0.7853981633974483i
So: math/base is okay for many ordinary complex values, excellent near 1+iε, but poor when the result has a very small nonzero real component. Chez’s complex algorithm is generally better for cancellation/overflow, except around the guarded branch point near 1.
While adding complex numbers to WebRacket I have found a potential problem
in the
math/baseimplementation ofatanh. When called with complexarguments, the accuracy isn't as good as Chez Scheme.
For the complex case the implementation of
atanhinmath/baseuses:This is mathematically correct, but not numerically stable.
The problem is how to improve the situation?
Option 1) improve the implementation in
math/base.Option 2) figure out how to call the Chez Scheme primitive
atanh.The section below was produced by Codex when I asked it examine the
accuracy of
atanhon complex numbers.Note: There might be similar problems with
asinhand andacosh./Jens Axel
Codex:
For
math/base, complex accuracy is not uniformly good. It uses the direct formula at [base-functions.rkt (line 104)](/Users/soegaard/Dropbox/GitHub/racket/racket/share/pkgs/math-lib/math/private/base/base-functions.rkt:104):That is mathematically right, but numerically fragile.
Measured against a 4096-bit bigfloat reference:
Examples:
The Chez source-level
cflatanhat 5_3.ss (line 761) is much more accuracy-oriented: it uses Kahan-style formulas and, whenlog1pis available through fllog1+ (line 412), it preserves those tiny real parts. In my model of that path, ordinary/tiny cases were typically 0-1 ulp.But Chez has its own weak spot near
z = 1+iε: it addsrho ≈ 2.98e-154to|imag|, so forε << rhoit clamps the behavior:So:
math/baseis okay for many ordinary complex values, excellent near1+iε, but poor when the result has a very small nonzero real component. Chez’s complex algorithm is generally better for cancellation/overflow, except around the guarded branch point near1.