Стабильный котангенс

Ответ 1

cot(x) = cos(x)/sin(x) должна быть более численно устойчивой вблизи π/2, чем cot(x) = 1/tan(x). Вы можете эффективно реализовать sincos на своих платформах.

Другая возможность - cot(x) = tan(M_PI_2 - x). Это должно быть быстрее, чем указано выше (даже если sincos доступно), но оно также может быть менее точным, так как M_PI_2 - это, конечно, только приближение трансцендентного числа π/2, поэтому разность M_PI_2 - x будет не быть точным для полной ширины a double мантиссы - на самом деле, если вам не повезло, у нее могут быть только несколько значимых бит.

Ответ 2

Если вы рассмотрите угол между двумя векторами (v и w), вы также можете получить котангенс следующим образом (используя Eigen :: Vector3d):

inline double cot(Eigen::Vector3d v, Eigen::Vector3d w) { 
    return( v.dot(w) / (v.cross(w).norm()) ); 
};

При тета-угле между v и w выше функция верна, потому что:

• | VxW | = | v |. | w |.sin (тета)
• против w = | v |. | w |.cos (тета)
• детская кроватка (тета) = cos (тета)/грех (тета) = (против w)/| vxw |

Ответ 3

TL; DR No.

Как правило, при поиске источников неточности следует в первую очередь заботиться о дополнениях и вычитаниях, которые могут привести к проблеме вычитания. Умножения и деления, как правило, безвредны для точности, кроме добавления дополнительной ошибки округления, но могут вызвать проблемы из-за переполнения и недостаточного значения при промежуточных вычислениях.

Ни один номер машины x может быть достаточно близко к кратным числам π/2, чтобы вызвать переполнение tan(x), поэтому tan(x) является четко определенным и конечным для всех кодировок с плавающей запятой для любой из чисел с плавающей запятой IEEE-754 форматы, и, соответственно, так же, как cot(x) = 1.0/tan(x).

Это легко продемонстрировать, выполнив исчерпывающий тест со всеми числовыми кодировками с float, поскольку исчерпывающий тест с использованием double невозможен, за исключением, возможно, самых больших суперкомпьютеров, существующих сегодня.

Используя математическую библиотеку с точной реализацией tan() с максимальной ошибкой ~ = 0,5 ulp, мы обнаруживаем, что при вычислении cot(x) = 1.0/tan(x) максимальная ошибка составляет менее 1,5 ulp, где дополнительная Ошибка по сравнению с самой tan() вносится ошибкой округления деления.

Повторяя этот исчерпывающий тест для всех значений с float с помощью cot(x) = cos(x)/sin(x), где sin() и cos() вычисляются с максимальной ошибкой ~ = 0,5 ulp, мы находим, что максимальная ошибка у cot() меньше 2,0 ulps, поэтому немного больше. Это легко объяснить наличием трех источников ошибок вместо двух в предыдущей формуле.

Наконец, cot(x) = tan (M_PI_2 - x) страдает от проблемы вычитания, упомянутой ранее, когда x около M_PI_2, и от проблемы, которая в арифметике с плавающей точкой конечной точности, M_PI_2 - x == M_PI_2 когда x достаточно мала по величине. Это может привести к очень большим ошибкам, из-за которых мы не получим правильных битов в результате.