Imagine a real-valued, discrete-time signal $h[n]$ with $N$ samples and a sampling period of $T_s = 1/f_s$. We know that the DFT of the signal $H[k]$ is conjugate symmetric because the signal is real-valued:

$H[N-k] = H^*[k], k=1,..,N-1$ and $H[0] \in R$

Now you want to apply a fractional delay $\tau [s]$ (i.e., any real-valued delay, not limited to an integer times the sampling period) in the frequency-domain:

$H_{\tau}[k] = e^{-\frac{j2 \pi k}{N}f_s\tau}H[k]$

Now verify symmetry again:

$H_{\tau}[N-k] = e^{-\frac{j2 \pi (N-k)}{N}f_s \tau}H[k] = e^{-j2 \pi f_s \tau}e^{+\frac{j2 \pi k}{N}f_s\tau}H[k]$

which means that the symmetry will only hold for $f_s \tau \in Z$ (i.e., the delay should be an integer times the sampling period). This also implies that the inverse DFT of my delayed spectrum will not be real-valued anymore, which actually makes no sense.

Therefore, my question is whether there's a straightforward way to apply a fractional delay in the frequency domain.

One way I can come up with is to apply the delay on one side of the spectrum, and then impose symmetry; however, this does not seem to yield the desired signal.

Another way is to multiply with the DFT of a fractionally delayed sinc-function, but then I am actually using a time-domain technique. Because the DFT of a sinc-function is actually a rectangular window, and delaying this sinc-function results in a phase shift in the frequency domain as described above, I understand that this way is some small correction to the first-described method which resulted in a complex-valued signal. However, I do not see what this small correction actually is in the frequency-domain, such that this technique does end up with a real-valued signal.

A third way I came up with is to do a DTFT instead of a DFT, apply the delay, and then sample this DTFT to end up with a DFT. However, I do not see how this works numerically.