sin(x)/x operator for LazyTensors

[haguettaz]

Adds sin(x)/x operator with the true value for the limit when x goes to 0. Warnings:

• the gradient is undefined in 0 EDIT: sinc'(0) = 0
• the pytorch version is not implemented, the comparison could fail if x == 0

Test Plan:

import pykeops
import torch
import math
from pykeops.torch import LazyTensor

device = 'cuda' 

x = torch.rand(1000, 1)*2*math.pi
y = x.data.clone()/math.pi # torch.sinc is defined as sin(pi*x)/pi*x
x = x.to(device)
y = y.to(device)
x.requires_grad = True
y.requires_grad = True

x_i = LazyTensor(x[:, None])
s1 = x_i.sinc().sum(0)
s2 = torch.sum(torch.sinc(y))
print("s1 - s2", torch.abs(s1 - s2).item())
assert torch.abs(s1 - s2) < 1e-3, torch.abs(s1 - s2)

s1.backward()
s2.backward()

print("grad_s1 - grad_s2", torch.max(torch.abs(x.grad - y.grad)).item())
assert torch.max(torch.abs(x.grad - y.grad)) < 1e-3

[bcharlier]

Hi @haguettaz ,

thank you for the PR. If I am correct, the derivative of sinc(x) = sin(x)/x could be extended to 0 when x==0:

[haguettaz]

Hi @bcharlier thanks for your comment. Indeed, I think you are right. The problem is I don't know how to define such a derivative (derivative by part) in the Sinc.h file. By the way, it appears pytorch implements this function in its new (unstable) version: torch.sinc, the unit test could be improved.

[bcharlier]

Hi @haguettaz ,

thanks again for your work! I suggest to use torch.sinc in the test.

[joanglaunes]

Hi and thanks again @haguettaz , A few remarks :

• In both PyTorch and NumPy, the convention for sinc is that it computes sin(pi x) / (pi x), so maybe we should do the same here ? (and maybe use your original sinxdivx for sin(x)/x ?). What do you think ?
• to properly implement the behavior of the gradient at x=0, I think the only way would be to implement a custom operation for the gradient, but I do not think it is worth the effort. In fact in PyTorch, I have checked and the sinc function returns NaN for gradient at 0, which is also the case with your implementation, so this is enough I think.

[haguettaz]

Hi @joanglaunes! You're welcome, it's a pleasure to help you.

1. Indeed, it is better to use the same convention as PyTorch and NumPy. Let's do that!
2. I agree with you, it is not necessary to implement the "true" gradient in this case. Maybe if PyTorch changes its implementation, we could think about it too. But at the moment, it seems to be enough.