Я пытаюсь реализовать arcsin в Python без использования какой-либо внешней библиотеки.
Вот мой код:
from time import process_time as pt
class TrigoCalc(metaclass=__readonly):
# This class evaluates various Trigonometric functions
# including Inverse Trigonometric functions
def __setattr__(self, name, value):
raise Exception("Value can't be changed")
@staticmethod
def asin(x):
'''Implementation from Taylor series
asin(x) => summation[(2k)! * x^(2k + 1) / (2^(2k) * (k!)^2 * (2k + 1))]
k = [0, inf)
x should be real
'''
# a0 = 1
# a1 = 1/(2*3)
# a2 = 1/2 * 3/(4*5)
# a3 = 1/2 * 3/4 * 5/(6*7)
# a4 = 1/2 * 3/4 * 5/6 * 7/(8*9)
# a5 = 1/2 * 3/4 * 5/6 * 7/8 * 9/(10*11)
# a6 = 1/2 * 3/4 * 5/6 * 7/8 * 9/10 * 11/(12*13)
# a7 = 1/2 * 3/4 * 5/6 * 7/8 * 9/10 * 11/12 * 13/(14*15)
# a8 = 1/2 * 3/4 * 5/6 * 7/8 * 9/10 * 11/12 * 13/14 * 15/(16*17)
# a9 = 1/2 * 3/4 * 5/6 * 7/8 * 9/10 * 11/12 * 13/14 * 15/16 * 17/(18*19)
# a10 = 1/2 * 3/4 * 5/6 * 7/8 * 9/10 * 11/12 * 13/14 * 15/16 * 17/18 * 19/(20*21)
# taking 10 coefficients for arriving at a common sequence
# N = n, D = n + 1; (N/D) --> Multiplication, number of times the coefficient number, n >= 1
start_time = pt()
coeff_list = []
NUM_ITER = 10000
for k in range(NUM_ITER):
if k == 0:
coeff_list.append(1)
else:
N = 1
D = N + 1
C = N/D
if k >= 2:
for i in range(k-1):
N += 2; D += 2
C = C * N/D
coeff_list.append(C)
__sum = 0
for k in range(NUM_ITER):
n = coeff_list[k] * math_utils.power(x, 2*k + 1) / (2*k + 1)
__sum += n
# Radian conversion to degrees
__sum = __sum/TrigoCalc.pi * 180
end_time = pt()
print(f'Execution time: {end_time - start_time} seconds')
return __sum
Результаты
Когда NUM_ITER составляет 60 (бесконечный ряд, повторяемый 60 раз), имеется значительная неточность в вычислениях на полюсе x = 1, тогда как x = 1/2 дает точность в 14 пунктов.
In [2]: TrigoCalc.asin(0.5)
Execution time: 0.0 seconds
Out[2]: 30.000000000000007
In [3]: TrigoCalc.asin(1)
Execution time: 0.0 seconds
Out[3]: 85.823908877692
Время выполнения в обоих прогонах незаметно .
Когда NUM_ITER равно 10000, то на x = 1 полюсе результат более точен, чем в предыдущем прогоне, но при x = 1/2 точность точно такая же.
In [4]: TrigoCalc.asin(0.5)
Execution time: 19.109375 seconds
Out[4]: 30.000000000000007
In [5]: TrigoCalc.asin(1)
Execution time: 19.109375 seconds
Out[5]: 89.67674183336727
Время выполнения в этих двух прогонах очень велико для этого типа вычислений.
Проблема
Как я могу сбалансировать код, чтобы он давал как минимум 1-балльную точность при x = 1 полюс в меньшем масштабе NUM_ITER?
Пожалуйста, не стесняйтесь давать предложения или обновления кода.
Python Ver: 3.7.7
РЕДАКТИРОВАТЬ: Изменения в коде для получения точных результатов с помощью ответа @ Joni
Заключение вычисления бесконечной серии в другую функцию внутри asin():
def asin(x):
def __arcsin_calc(x):
# ....
# Computations
# ....
# Removing the radian to degree conversion from this function
return __sum
Добавление ограничений к x с помощью новой функции внутри asin() для избежать медленной сходимости:
if -1.0 <= x < -0.5:
return -(TrigoCalc.pi/2 - __arcsin_calc(math_utils.power((1 - x*x), 0.5))) / TrigoCalc.pi * 180 # Radian to Degree conversion
elif -0.5 <= x <= 0.5:
return __arcsin_calc(x)/TrigoCalc.pi * 180
elif 0.5 < x <= 1.0:
return (TrigoCalc.pi/2 - __arcsin_calc(math_utils.power((1 - x*x), 0.5))) / TrigoCalc.pi * 180
else:
raise ValueError("x should be in range of [-1, 1]")
Результаты:
In [2]: TrigoCalc.asin(0.99)
Execution time: 0.0 seconds
Out[2]: 81.89022502527023
In [3]: math.asin(0.99)/TrigoCalc.pi*180
Out[3]: 81.89038554400582
In [4]: TrigoCalc.asin(1)
Execution time: 0.0 seconds
Out[4]: 90.0
In [5]: math.asin(1)/TrigoCalc.pi*180
Out[5]: 90.0