Ответ
np.sin
в общем случае будет максимально точным, учитывая точность переменных double
(то есть 64-битных float
), в которых хранятся входные, выходные и промежуточные значения. Вы можете получить разумную меру точности np.sin
, сравнив ее с версией произвольной точности sin
из mpmath
:
import matplotlib.pyplot as plt
import mpmath
from mpmath import mp
# set mpmath to an extremely high precision
mp.dps = 100
x = np.linspace(-np.pi, np.pi, num=int(1e3))
# numpy sine values
y = np.sin(x)
# extremely high precision sine values
realy = np.array([mpmath.sin(a) for a in x])
# the end results are arrays of arbitrary precision mpf values (ie abserr.dtype=='O')
diff = realy - y
abserr = np.abs(diff)
relerr = np.abs(diff/realy)
plt.plot(x, abserr, lw=.5, label='Absolute error')
plt.plot(x, relerr, lw=.5, label='Relative error')
plt.axhline(2e-16, c='k', ls='--', lw=.5, label=r'$2 \cdot 10^{-16}$')
plt.yscale('log')
plt.xlim(-np.pi, np.pi)
plt.ylim(1e-20, 1e-15)
plt.xlabel('x')
plt.ylabel('Error in np.sin(x)')
plt.legend()
Выход:
Таким образом, разумно сказать, что как относительные, так и абсолютные ошибки np.sin
имеют верхнюю границу 2e-16
.
Лучший ответ
Существует отличный шанс, что если вы сделаете increment
достаточно маленьким, чтобы ваш подход был точным, ваш алгоритм будет слишком медленным для практического использования. Подходы к решению стандартных уравнений не будут работать для вас, поскольку у вас нет стандартной функции. Вместо этого у вас есть неявная многозначная функция. Вот пример общего подхода для получения всех решений этого вида уравнения:
import matplotlib.pyplot as plt
import numpy as np
import scipy.optimize as spo
eps = 1e-4
def func(x, a):
return a*np.sin(x) - x
def uniqueflt(arr):
b = arr.copy()
b.sort()
d = np.append(True, np.diff(b))
return b[d>eps]
initial_guess = np.arange(-9, 9) + eps
# uniqueflt removes any repeated roots
roots = uniqueflt(spo.fsolve(func, initial_guess, args=(10,)))
# roots is an array with the 7 unique roots of 10*np.sin(x) - x == 0:
# array([-8.42320394e+00, -7.06817437e+00, -2.85234190e+00, -8.13413225e-09,
# 2.85234189e+00, 7.06817436e+00, 8.42320394e+00])
x = np.linspace(-20, 20, num=int(1e3))
plt.plot(x, x, label=r'$y = x$')
plt.plot(x, 10*np.sin(x), label=r'$y = 10 \cdot sin(x)$')
plt.plot(roots, 10*np.sin(roots), '.', c='k', ms=7, label='Solutions')
plt.ylim(-10.5, 20)
plt.gca().set_aspect('equal', adjustable='box')
plt.legend()
Выход:
Вам нужно настроить initial_guess
в зависимости от вашего значения a
. initial_guess
должно быть не меньше фактического количества решений.