Najlepsza odpowiedź
Masz na myśli funkcję obliczającą \ frac {\ sin x} {x} dla x \ ne 0 i 1 dla x = 0?
Oczywiście łatwo byłoby po prostu zaprogramować coś za pomocą instrukcji takiej jak
function result = mySinc(x)
\% My implementation of the sinc function.
if x == 0
result = 1;
else
result = sin(x)/x;
end
end
Ale w rzeczywistości może to być mniej niż idealny. Jasne, unika się osobliwości i ma prawidłową strukturę matematyczną, ale oczywiście musimy pamiętać, że wartości w x nie są liczbami „t”: to są wzorce bitów reprezentowane przez przełączniki w sprzęcie i tylko przybliżone liczby. Działają całkiem dobrze przez większość czasu, ale zdarzają się sytuacje, w których coś po prostu nie działa. Na przykład podczas oceniania a + log(1+x) dla małych x lub podczas oceny a-b, gdy a i b mają bardzo podobne wartości. może stracić dużo precyzji, jeśli ocenisz to drugie naiwnie (chociaż czasami możesz nie mieć wyboru, niestety). Dzielenie przez małe liczby również może powodować problemy, chociaż na szczęście oba \ sin x i x są porównywalne, gdy | x | \ ll 1. Może nie byłby to problem, ale w przypadku arytmetyki zmiennoprzecinkowej rzadko jest dobrą praktyką posiadanie kodu takiego jak if x == foo.
Co zrobić w tym przypadku? Rozważmy szereg Taylora \ frac {\ sin x} {x} wyśrodkowany na x = 0. Nie zapiszę tego, ale na pewno wiesz, jak to obliczyć.
Następnie zamień test if x == 0 na if abs(x) < sqrt(eps) , gdzie w Matlabie eps zwraca maszynę epsilon (porównywalną do 10 ^ {- 16} w większości dzisiejszych systemów) i zastępuje wiersz result = 1; by result = ...;, gdzie ... będzie obliczeniem szeregu Taylora obciętego po odpowiednim porządku (zwykle około 2 powinno ale powinieneś sam sprawdzić skalę najwyższego terminu zamówienia, które wybrałeś dla x \ sim 10 ^ {- 8} i czy wynik będzie tak dokładny, jak potrzebujesz, czy nie).
fajną rzeczą jest to, że to, co stracisz w szybkości wykonywania tego kodu, prawdopodobnie zwiększysz dokładność, szczególnie podczas obliczania całki sinusoidalnej \ operatorname {Si} (x) = \ int\_0 ^ xf (t) dt, gdzie f (t ) jest ciągłym rozszerzeniem funkcji \ operatorname {sinc} na wszystkie rzeczywiste jako jej domenę.