[Numpy-discussion] Converting np.sinc into a ufunc

[Joshua Wilson]

Re Ralf's question:

> Can you quantify the precision improvement (approximately)?

On one level you'll get a large decrease in relative error around the
zeros of the sinc function because argument reduction is being done by
a number which is exactly representable in double precision (i.e. the
number 2) versus an irrational number. For example, consider:

>>> import numpy as np
>>> import mpmath
>>> x = 1 + 1e-8
>>> (np.sinc(x) - mpmath.sincpi(x))/mpmath.sincpi(x)

On master that will give you

mpf('-5.1753363184721223e-9')

versus

mpf('1.6543612517040003e-16')

on the new branch.

*But* there are some caveats to that answer: since we are close to the
zeros of the sinc function the condition number is large, so real
world data that has already been rounded before even calling sinc will
incur the same mathematically unavoidable loss of precision.

On Wed, May 22, 2019 at 6:24 PM Charles R Harris
<charlesr.harris at gmail.com> wrote:
>
>
>
> On Wed, May 22, 2019 at 7:14 PM Marten van Kerkwijk <m.h.vankerkwijk at gmail.com> wrote:
>>
>> On a more general note, if we change to a ufunc, it will get us stuck with sinc being the normalized version, where the units of the input have to be in the half-cycles preferred by signal-processing people rather than the radians preferred by mathematicians.
>>
>> In this respect, note that there is an outstanding issue about whether to allow one to choose between the two: https://github.com/numpy/numpy/issues/13457 (which itself was raised following an inconclusive PR that tried to add a keyword argument for it).
>>
>> Adding a keyword argument is much easier for a general function than for a ufunc.
>>
>
> I'd be tempted to have two sinc functions with the different normalizations. Of course, one could say the same about trig functions in both radians and degrees. If I had to pick one, I'd choose sinc in radians, but I think that ship has sailed.
>
> Chuck
> _______________________________________________
> NumPy-Discussion mailing list
> NumPy-Discussion at python.org
> https://mail.python.org/mailman/listinfo/numpy-discussion