Suppose you are given two prime numbers, $p<q$ such that $p\mid\varphi(q)$.
My question can be stated in two equivalent manners:
Let $G$ be the group of units $\mathbb{Z}_q^\times. $Is there an efficient algorithm that will enable me to find an element $g\in G$ (other than $1$) such that $g^p\equiv1$?
Is there an efficient algorithm that will let me solve the cyclotomic polynomial $\Phi_p(g)$ over the finite field $\mathbb{F}_q$?
Obviously, brute force can be used to find solutions, so my question is about an algorithm that can solve these questions in a manner that is significantly faster than brute force.
For $p=2$, the problem is easy, and the solution is that $g=q-1$. The problem is more difficult for other values of $p$.
Let $H$ be the cyclic group $\mathbb{Z}_{q-1}$. Finding an element $h_0\in H$ of order $p$ is trivial, as $h_0=\frac{q-1}{p}$ is a solution. However, even though $G\cong H$, this doesn't help us unless we can also find an isomorphic mapping $\psi: H\rightarrow G$. However, if we randomly guess possible values $g_0\in G$, we can stop as soon as we find an element of either order $p$ (giving us an immediate solution of $g=g_0$), or of order $q-1$ (giving us a mapping $\psi$, and a resulting solution $g=g_0^{h_0}$). As there are $\varphi(p)$ elements in $G$ of order $p$ and $\varphi(q-1)$ elements of $G$ of order $q-1$, this guess and check algorithm, on average, will yield a solution in $\frac{\varphi(q)-2}{\varphi(p)+\varphi(q-1)}$ guesses (the $-2$ in the numerator comes from the fact that we don't need to bother checking the elements $1,q-1$, as we know their orders are 1 and 2, respectively). For any fixed value of $p$, this grows approximately linearly with $q$.
For $p=3$ or $p=5$, we can use the quadratic formula to solve for $\Phi_3$ or the quartic formula to solve for $\Phi_5$ respectively. This gives us that for $p=3$, we have an example solution of $g=\frac{-1+\sqrt{-3}}{2}$ and for $p=5$, we have an example solution of $y=\sqrt{5}, g=\frac{y+1}{4}+\sqrt{\frac{5-y}{8}}$. There are known algorithms for solving square roots in $\mathbb{F}_q$ that grow in runtime approximately linearly with the number digits of $q$, which makes them significantly faster for large $q$.
However, I can't find a way to extend this methodology to other values of $p$. In order for the method to be applicable, we require that $\Phi_p$ have a known solution with radicals, and that every type of radical used in that explicit solution has an algorithm (more efficient than brute force) to solve over a finite field. From a quick Google search, I failed to find any such algorithms for solving cube roots, and I also failed to find an explicit solution in radicals for $\Phi_7$ (although I did find articles stating that the Galois groups of the splitting fields of cyclotomic polynomials are all abelian, which means that such explicit solutions must exist).
EDIT: I found an article giving the explicit solution for $\Phi_7$ in radicals, but didn't explain how it got the solution, other than by saying that the methodology does not extend to higher primes. The solution depends on solving four radicals: $x=\sqrt{-3}, y=\sqrt[3]{14x}, z=\sqrt[3]{-7x}, w=\sqrt{(\frac{y+z-1}{3})^2-4}$. The explicit solution can then be written as $\frac{y+z+3w-1}{6}$.