[Bug libfortran/93871] COTAN is slow for complex types

[thenlich at gcc dot gnu.org]

https://gcc.gnu.org/bugzilla/show_bug.cgi?id=93871

--- Comment #19 from Thomas Henlich <thenlich at gcc dot gnu.org> ---
Regarding the following:

https://stackoverflow.com/questions/3738384/stable-cotangent#56864824
====
Is there a more stable implementation for the cotangent function than return
1.0/tan(x)?

No.

No machine number x can get close enough to multiples of π/2 to cause tan(x) to
overflow, therefore tan(x) is well-defined and finite for all floating-point
encodings for any of the IEEE-754 floating-point formats, and by extension, so
is cot(x) = 1.0 / tan(x).

...

Using a math library with an accurate implementation of tan() with a maximum
error of ~= 0.5 ulp, we find that computing cot(x) = 1.0 / tan(x) incurs a
maximum error of less than 1.5 ulp, where the additional error compared to
tan() itself is contributed by the rounding error of the division.

Repeating this exhaustive test over all float values with cot(x) = cos(x) /
sin(x), where sin() and cos() are computed with a maximum error of ~= 0.5 ulp,
we find that the maximum error in cot() is less than 2.0 ulps, so slightly
larger. This is easily explained by having three sources of error instead of
two in the previous formula.
====

it may be best (and fastest) to just implement it like above: cot=1/tan