Density of a set of numbers dividing a fixed number with polynomial exponent

Question

Fix a positive integer $a>1$ and let $f\in\mathbb{Z}[x]$ be a polynomial with positive leading coefficient. We define a set $S$ of positive integers, $$ S=\{n\in\mathbb{Z}^+:n\mid a^{f(n)}-1\}. $$ Question is, can we compute the density of $S$, that is, $$ \lim_{n\to\infty}\frac{|S\cap \{1,2,\dots,n\}|}{n}? $$ I believe that the answer should be $0$, yet I don't have a proof.

Answer 1

Fix $\varepsilon>0$ and choose large prime $q>5\varepsilon^{-1} \cdot \deg(f)$ such that $f$ is non-trivial modulo $q$ (then $f$ has at most $\deg(f)$ roots modulo $q$) and additionally such that $a$ is not a $q$-th perfect power. Say that a prime $p>q$ is $q$-appropriate if $a$ is not $q$-th power modulo $p$. By Chebotarev density theorem the $q$-appropriate primes have positive density. Thus we may choose $q$-appropriate primes $p_1,\dots,p_m$ such that $\prod_{i=1}^m (1-1/p_i)<\varepsilon/5$. Note that if $n$ divides $a^{f(n)}-1$, then either $q$ divides $f(n)$ (density of such $n$ does not exceed $\deg(f)/q<\varepsilon/5$), or neither $p_i$ divides $n$. [Indeed, if $p_i|n|a^{f(n)}-1$, but $q$ does not divide $f(n)$, we may find $s$ such that $sf(n)+1$ is divisible by $q$, that yields $a\equiv a^{1+sf(n)} \pmod {p_i}$, so $a$ is a $q$-th power modulo $p_i$, a contradiction.] The density of such $n$'s is exactly $\prod_{i=1}^m (1-1/p_i)<\varepsilon/5$. So $S$ is covered by finitely many arithmetic progressions with sum of densities less than $\varepsilon$.