Thank you for this great work. To the extent I implemented some of the missing trig functions in Complex, I followed the implementations in Complex.js , which seemed well worked out and supported. However I am sure Boost would be much better engineered still. So +1
On Tue, Dec 10, 2019 at 7:42 AM Alex Herbert <alex.d.herb...@gmail.com> wrote: > The latest round of fixes on Complex have made it fully C99 compliant > for a range of -5 to +5 in the real and imaginary components. I have yet > to finish the testing of overflow/underflow conditions. Here I report on > the standard range data. > > I have updated the test of Complex to load test data from file > resources. These files have been produced using GNU gcc. But the data > can be swapped for another reference source. > > Currently the test is using very high tolerances specified in units of > least precision (ULPs). Dropping in another file using low precision > data (not the full 17 fractional digits of a double) would break the > test. So this may have to be revised in the future if other references > are to be used. > > Measuring the difference between Complex and the reference data with the > units of least precision (ULP) delta between them shows: > > function mean +/- sd (n=count) [min, max] Median=x > > acos 3.94 +/- 5.03 (n=484) [ 0, 36] Median= 2 > acosh 3.94 +/- 5.03 (n=484) [ 0, 36] Median= 2 > asinh 0.05 +/- 0.23 (n=242) [ 0, 1] Median= 0 > atanh 0.84 +/- 2.14 (n=242) [ 0, 17] Median= 0 > cosh 0.09 +/- 0.31 (n=242) [ 0, 2] Median= 0 > sinh 0.08 +/- 0.30 (n=242) [ 0, 2] Median= 0 > tanh 0.82 +/- 2.30 (n=242) [ 0, 34] Median= 0 > exp 0.06 +/- 0.31 (n=484) [ 0, 2] Median= 0 > log 0.10 +/- 0.37 (n=484) [ 0, 3] Median= 0 > sqrt 0.02 +/- 0.16 (n=484) [ 0, 1] Median= 0 > multiply 0.00 +/- 0.00 (n= 32) [ 0, 0] Median= 0 > divide 1.03 +/- 1.75 (n= 32) [ 0, 7] Median= 1 > pow 1.12 +/- 2.61 (n= 32) [ 0, 9] Median= 0 > > If we only include those with a delta above 1 ULP (i.e. so that there is > at least one number between the two): > > acos 6.54 +/- 5.38 (n=274) [ 2, 36] Median= 5 > acosh 6.54 +/- 5.38 (n=274) [ 2, 36] Median= 5 > asinh NaN +/- NaN (n= 0) [NaN, NaN] Median=NaN > atanh 4.56 +/- 3.96 (n= 34) [ 2, 17] Median= 3 > cosh 2.00 +/- 0.00 (n= 2) [ 2, 2] Median= 2 > sinh 2.00 +/- 0.00 (n= 2) [ 2, 2] Median= 2 > tanh 3.22 +/- 5.30 (n= 36) [ 2, 34] Median= 2 > exp 2.00 +/- 0.00 (n= 9) [ 2, 2] Median= 2 > log 3.00 +/- 0.00 (n= 4) [ 3, 3] Median= 3 > sqrt NaN +/- NaN (n= 0) [NaN, NaN] Median=NaN > multiply NaN +/- NaN (n= 0) [NaN, NaN] Median=NaN > divide 5.00 +/- 2.31 (n= 4) [ 3, 7] Median= 5 > pow 5.83 +/- 3.06 (n= 6) [ 2, 9] Median= 7 > > This mainly shows a count of real differences ignoring floating-point > round-off. > > These show that Complex is doing quite well for all but: > > acos > acosh > atanh > tanh > > acosh is implemented using a trigonomic identity with acos and the two > are equally bad. Fixing acos would fix this too. > > I am not sure what can be done for tanh. This uses the formula: > > tan(a + b i) = sinh(2a)/(cosh(2a)+cos(2b)) + i [sin(2b)/(cosh(2a)+cos(2b))] > > It is implemented entirely using the Math library. A histogram of the > failures show that there is mainly one anomaly with a ULP delta of 34: > > tanh 2 23 > tanh 3 12 > tanh 34 1 > > For now this can be left for further investigation. More data may show > that this error is an outlier. > > The data is: (0.0,1.5).tanh(). Perhaps there is a trigonomic identity to > use when the real (or imaginary) component is zero. > > So how to fix acos and atanh? > > atanh uses divide and log on a Complex. > > acos uses multiply, sqrt and log on a Complex. > > Thus the differences between Complex and the standard are in those > methods that are using the Complex object to perform part of the > computation. > > I note that asin uses asinh which uses multiply, sqrt and log on a > Complex and this does not perform badly. The method is very similar to > acos. This may be due to the range of the test data or the > implementation using only positive component parts to preserve the > conjugate and odd function equalities. > > The C++ boost library has implementations for asin, acos and atanh. > These have overflow and underflow protection and an efficient > computation in a 'normal' value range. I will investigate using those > methods to see if they make a difference to the error. > > I would like to get this error down to under 1 ULP on average with a > lower maximum. > > Then I will look at the boundary cases for finite numbers which have > overflow or underflow during parts of the computation. I am sure that > the functions using Math.sinh and Math.cosh which asymptote to infinity > require some overflow protection. These are: > > cosh > sinh > tanh > > Alex > > >
