Вычисление арксинуса с высокой точностью Опции

[snav]

Ребят, возникла такая задача на работе: нужно научиться вычислять тригонометрические функции с высокой точностью. В стандартном 8-байтном типе double слишком мало разрядов, поэтому библиотечные функции не годятся. Хочу использовать 16-байтный тип decimal в .NET. Функции Sin, Cos, Tan, Atan уже написал, используя разложение в степенные ряды. Осталось придумать, как запрограммировать арксинус и арккосинус.

Стандартный степенной ряд arcsin(x) = x + 1/6 * x^3 + 3/40 * x^5 +...
сходится более менее быстро только при |x| < 1/2. А при |x| > 1/2 скорость сходимости резко падает. Сейчас вычисляю арксинус через численное решение уравнения sin(x)=a, но это очень медленно.
Чтобы придумать для быстрого вычисления в диапазоне |x| > 1/2?

Может быть, выполнить какое-то математическое преобразование, чтобы от вычисления арксинуса в области |x| > 1/2 перейти к вычислению арксинуса в области |x| < 1/2. Но преобразование не должно приводить к существенной потере точности.

[kikotpb]

Пусть 1/2<=x<=1, arcsin x=a, П/6 <=a<=П/2. Положим b=a-П/3, |b|<=П/6.
sin b= sin(a-П/3)=sin a cos П/3-cos a sin П/3 =1/2 sin a - sqrt(3)/2 cos a = x/2-(sqrt(1-x^2))*sqrt(3)/2;
b=arcsin (x/2-(sqrt(1-x^2))*sqrt(3)/2). Т.к. |b|<=П/6, то |x/2-(sqrt(1-x^2))*sqrt(3)/2|<=1/2, а это Вы уже умеете. Наконец, a=arcsin x=b+П/3.
Правда, приходится вычислять корни, но тут я не знаю, насколько оно критично.

[snav]

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

[0]

Ну вообще по классике замечательно считается арктангенс а через него арксинусы и арккосинусы легко считаются.

[snav]

Пусть 1/2<=x<=1, arcsin x=a, П/6 <=a<=П/2. Положим b=a-П/3, |b|<=П/6.
sin b= sin(a-П/3)=sin a cos П/3-cos a sin П/3 =1/2 sin a - sqrt(3)/2 cos a = x/2-(sqrt(1-x^2))*sqrt(3)/2;
b=arcsin (x/2-(sqrt(1-x^2))*sqrt(3)/2). Т.к. |b|<=П/6, то |x/2-(sqrt(1-x^2))*sqrt(3)/2|<=1/2, а это Вы уже умеете. Наконец, a=arcsin x=b+П/3.

Попробовал реализовать эту идею на практике, вроде работает.
Но я сделал немного по-другому. Пусть y=arcsin(x). Сначала получаем грубую оценку значения арксинуса y0 с помощью обычной библиотечной функции для типа double, а затем уточняем это значение, вводя вспомогательную переменную z = y - y0.
В упрощенном виде (без проверок и обработки узких мест) это выглядит так:

CODE

decimal Asin(decimal x)
{
    decimal y0 = (decimal)Math.Asin((double)x);

    decimal sin_y  = x;
    decimal cos_y  = Sqrt(1-x*x);
    decimal sin_y0 = Sin(y0);
    decimal cos_y0 = Cos(y0);
    decimal sin_z  = sin_y * cos_y0 - cos_y * sin_y0;

    decimal z = _arcsin(sin_z);   // Вычисление арксинуса через степенной ряд, при этом на практике z = sin_z
    decimal y = z + y0;

    return y;
}

Удобством этого способа является то, что фактически отпадает необходимость вычислять arcsin(sin_z), так как z очень мало (~1e-15) и в этом случае z=sin_z с точностью до величины порядка z^2 = 1e-30, т.е. ошибка меньше, чем разрядность числа decimal. На всякий случай я всё равно вычисляю z через ряд, но метод _arcsin(sin_z) сразу обнаруживает обнуление второго члена ряда и возвращает sin_z.

К недостаткам можно отнести необходимость вычисления Cos(y0), Sin(y0) и Sqrt(1-x*x), что делает расчет арксинуса более долгим по сравнению с расчетом синуса, но не значительно. Корень с помощью итерационной процедуры Ньютона находится на удивление быстро.

Остается, правда, открытым вопрос о точности такого способа вычисления, но я пока не обнаружил каких-то больших ошибок. При проверке на множестве чисел в диапазоне -1...1, расхождение между x и Sin(Asin(x)) не превышало 1e-27, т.е. ошибка в пределах одного десятичного разряда числа decimal, что можно считать неплохим результатом, учитывая что с такой точностью я вычислял значение корня.

[nik_vic]

Ну вообще по классике замечательно считается арктангенс а через него арксинусы и арккосинусы легко считаются.

http://algorithm.narod.ru/el/arctan.html