Peter asked here "Can a number be equal to the sum of the squares of its prime divisors?" and, it seems clear that if $$n=p_1^{a_1}\cdots p_k^{a_k},$$ and $$f(n):=p_1^2+\cdots+p_k^2$$ that then $n=f(n)$ only if $n$ is a square of a prime.

I almost proved that, but did not yet find so good and simple bounds to show that indeed it is possible only if $n$ is a square of a prime.

But, despite the fact I still do not have a proof, I decided to ask a little more general question:

If $$n=p_1^{a_1}\cdots p_k^{a_k},$$ is such that $$n=b_1p_1^2+...+b_kp_k^2$$ for some $b_1,...b_k \in \mathbb N$, then what is the density of the set of all such $n$ in $\mathbb N$?

I am also interested in how often do you find more than the one set $(b_1,...b_k)$ that solves $p_1^{a_1}...p_k^{a_k}=b_1p_1^2+...+b_kp_k^2$? That is, how often a function $n \to (b_1,...b_k)$ is multi-valued?

What is the biggest $n$ for which you did not find a solution?

Is the sequence of $n$´s for which there is no solution in OEIS?