Roots of unity in non-Abelian groups: when do they form subgroups?

I am going to use the notation $\sqrt[n]{G}$ instead. In this answer all groups will be finite.

If $G$ has a Sylow $p$-subgroup $P$ which is not normal, then the elements of order dividing $|P|$ can't form a subgroup (since they consist precisely of the conjugates of $P$). Hence in order for $\sqrt[n]{G}$ to always be a subgroup, every Sylow subgroup needs to be normal.

Lemma: Let $G$ be a finite group all of whose Sylow subgroups are normal. Then $G$ is the direct product of its Sylow subgroups.

Proof. Let $g$ be an element of a Sylow $p$-subgroup and $h$ be an element of a Sylow $q$-subgroup, $p \neq q$. Then $ghg^{-1} h^{-1}$ lies in a Sylow $q$-subgroup (since $ghg^{-1}$ is in a Sylow $q$-subgroup) but also in a Sylow $p$-subgroup (since $hg^{-1} h^{-1}$ lies in a Sylow $p$-subgroup), and the intersection of any two such subgroups is trivial. Hence $gh = hg$. It follows that the obvious map from the direct product of the Sylow subgroups of $G$ to $G$ is an injective homomorphism between two groups of the same size, hence an isomorphism. $\Box$

Write $G = P_1 \times ... \times P_k$ where $P_i$ is the Sylow $p_i$-subgroup ($p_i$ the primes in increasing order) and let $n = \prod p_i^{e_i}$. Then

$$\sqrt[n]{G} = \prod \sqrt[p_i^{e_i}]{P_i}.$$

Consequently, $\sqrt[n]{G}$ is always a subgroup of $G$ if and only if $\sqrt[p_i^{e_i}]{P_i}$ is always a subgroup of $P_i$.

So the problem reduces to the case of non-abelian $p$-groups, and here I have no idea what can happen in general, but here are infinitely many examples for every $p \ge 3$: the discrete Heisenberg group $H_3(\mathbb{F}_{p^n})$ is non-abelian of order $p^{3n}$, and every element can be written in the form $I + N$ where $N$ is a $3 \times 3$ nilpotent matrix over $\mathbb{F}_{p^n}$. Hence

$$(I + N)^p = I + p N + {p \choose 2} N^2 = I$$

and we conclude that every element has order dividing $p$.