[emrys53]

Fixes #1432.
Supersedes #1434, #1435, #1436, #1437, #1438, #1439, #1440.

The implementations are mostly taken from FreeBSD implementation (scalar, converted to simd version) and some from Julia. Main idea is that we always try to sum up variables that have the same sign to avoid cancellations and in some cases we check for some special fixed bounds so that the term won't become Inf. I also added tests, that tests for the inputs mentioned in #1432.

Performance: ns per element, master vs this PR

Lower is faster; the last column of each architecture is PR / master (below 1 means the PR is faster).

function	type	SSE2 master	SSE2 PR	SSE2 PR/master	AVX2+FMA master	AVX2+FMA PR	AVX2+FMA PR/master	AVX512F master	AVX512F PR	AVX512F PR/master
asin	float	11.33	8.53	0.76	5.52	4.33	0.78	3.11	2.82	0.90
asin	double	28.99	25.79	0.90	14.84	13.33	0.90	12.78	11.46	0.90
acos	float	11.25	8.51	0.76	5.54	4.35	0.79	3.14	2.76	0.87
acos	double	28.95	26.20	0.90	14.81	13.36	0.91	12.67	11.54	0.91
asinh	float	11.26	8.36	0.74	5.60	4.35	0.78	3.12	2.75	0.88
asinh	double	29.13	25.71	0.90	14.83	13.42	0.91	12.76	11.49	0.90
acosh	float	11.18	8.47	0.78	5.54	4.25	0.77	3.14	2.79	0.88
acosh	double	28.81	26.12	0.90	14.81	13.43	0.91	12.79	11.53	0.90
tan	float	14.99	11.51	0.76	9.11	6.45	0.71	9.81	6.29	0.64
tan	double	23.41	20.39	0.88	7.76	8.00	1.04	6.23	6.64	1.07
tanh	float	14.46	10.97	0.78	9.26	6.30	0.68	10.14	6.49	0.63
tanh	double	22.87	20.18	0.87	7.47	8.11	1.07	5.70	6.80	1.20
log1p	float	5.08	5.42	1.07	2.39	2.67	1.05	1.48	1.49	1.04
log1p	double	13.81	14.51	1.05	7.00	7.48	1.04	5.06	5.13	1.01
atan2	float	1.75	1.76	1.00	0.77	0.81	1.05	0.55	0.57	1.04
atan2	double	4.55	4.98	1.10	2.31	2.24	1.03	2.35	2.30	0.98

• Machine: Intel Xeon w5-3435X. Compiler: GCC 14.2, -O3 -DNDEBUG.

[emrys53]

All failing CI jobs except one use libc++, whose trigonometric functions use different method and are less accurate than the current xsimd for large |z|. I checked some of the values with "exact" results and xsimd is a lot closer to "exact" results than what libc++ provides. I don't know how to fix this issue. Maybe instead of testing against std:: implementations, we can test it against precomputed values for each case so it doesn't change between different environments.

The one exception is clang-cl x64 /fp:fast. The test hits log1p(-1), which should be -inf. With fast-math on, xsimd's log doesn't handle infinities, so it returns a finite value instead. I think with -ffast-math we can assume there should be no NaN/Inf for arguments/results so in my opinion xsimd's approach seems better than what clang does. Maybe we can disable these -inf cases for this specific envrionment.

[DiamonDinoia]

We can change the tests to be accurate to machine precision multiplied by the conditioning of the problem.
I would say that we need to be backward stable (threfeten & bau). If we account for the conditioning of the problem we lost roughly one digit for every power of 10 of the magnitude of the argument.

Once we do that we can just measure the relative error point-wise to that accuracy.

fast-math is best effort user that enable it know that they cannot distinguish inf from a large value.

[DiamonDinoia]

I do not like testing against generated data from python.
Nor I like a 3rd party dependency for those.
How hard is to hardcode the corner cases we care? There should be very few of those.

[DiamonDinoia]

Also, can you check the license of the code you took? I do not want to violate any of them. We can just implement the textbook algorithm for these things and fit the coefficients with Remez.

[DiamonDinoia]

Over engineered tests :)

[emrys53]

Also, can you check the license of the code you took? I do not want to violate any of them. We can just implement the textbook algorithm for these things and fit the coefficients with Remez.

Licences for FreeBSD are BSD-2 and for Julia it is MIT. I will add the Licences and credit on top of each function in the next commit.

[DiamonDinoia]

since there are citations I do not need to restate the equations. Also I can read the equations from the code itself.

[emrys53]

since there are citations I do not need to restate the equations. Also I can read the equations from the code itself.

I tried to keep the comments to minimum and only explain math and why some steps are necessary and why is it done this way to make it easier to review (and also myself to keep track of what is going on). I tried to keep them to minimum but if you think they are too much i can remove them. But in my opinion it makes the code more clear and it is not too much comment.