Because a coefficient of $ x^{2d}$ in a polynomial $ (x+\frac{1}{x})^{2k}$ is $ {2k}\choose {k+d}$, We should calculate a coefficient of $ x^{2d}$ in $ P(x)=\sum_{k=0}^{p-1}(x+\frac{1}{x})^{k} = \frac{(x+\frac{1}{x})^{2p}-1}{(x+\frac{1}{x})^{2}-1}=\frac{(x^{2}+1)^{2p}}{x^{2p+1}-x^{2p}+x^{2p-1}}-\frac{x}{x^{2}-x+1}$. But each coefficient of $ Q(x)^{p}$ is congrent to each coeffient of $ Q(x^{p}) \pmod p$ when $ Q$ is a polynomial since $ p\choose 1$ $ \sim $ $ p\choose {p-1}$$ \equiv 0 \pmod p$. Therefore $ P(x)\equiv_{p} \frac{(x^{2p}+1)^{2}}{x^{2p+1}-x^{2p}+x^{2p-1}}-\frac{x}{x^{2}-x+1}=\left( \frac{x}{x^{2}-x+1}\right)(x^{2p}+1+\frac{1}{x^{2p}})$ Let $ Q(x)=\sum_{i=0}^{\infty}a_{i}x^{i}$ when $ a_{3i}=0$, $ a_{6i+1}=a_{6i+2}=1$, $ a_{6i+4}=a_{6i+5}=-1$ for integer $ i$. Then $ a_{n+2}=a_{n+1}-a_{n}$. So following is holds : $ Q(x)-x=xQ(x)-x^{2}Q(x)$, which we get $ Q(x)=\frac{x}{x^{2}-x+1}$. Now we know that $ f(d)=[x^{2d}]P(x)= [x^{2d}]Q(x)+[x^{2d+2p}]Q(x) = a_{2d}+a_{2d+2p}$. Remains to show is that $ f(d)=r$, it can be proved by just case-by-case calculation.