When is the sum of divisors a perfect square? For $n=3$, $\sigma(n)=4$, a perfect square. Calculating further was not yielding positive results.
I was wondering is there a way to find all such an $n$, like some algorithm?
We know that if $n=p_1^{k_1}p_2^{k_2}\dots p_r^{k_r}$ then $\sigma(n)=\frac {p_1^{k_1+1}-1}{p_1-1} \frac{p_2^{k_2+1}-1}{p_2-1}\dots \frac {p_r^{k_r+1}-1}{p_r-1}$. So we have to find $n$ such that $\sigma(n)=\frac {p_1^{k_1+1}-1}{p_1-1} \frac{p_2^{k_2+1}-1}{p_2-1}\dots \frac {p_r^{k_r+1}-1}{p_r-1}=t^2$ for some integer $t$. If each $\frac{p_i^{k_i+1}-1}{p_i-1}$ is a square we will get some such numbers (not all though), but let us consider $\frac{p_i^{k_i+1}-1}{p_i-1}=m^2$. 
Working a example, If I take $p_i=17 $ , (or any prime such that $p_i-1$ is a perfect square), then if I can find a $k_i$ such that ${p_i^{k_i+1}-1}=16.m^2$ I will get an $n$, so now when can this happen, is something I have to try out. 
But out of curiosity I am asking has this been worked out before. Problem looks broad to me, Is there a solution to this problem?
 A: I came up with an algorithm, decades ago. Kap was interested in solving $\sigma(x^3) = y^2,$ where the next simplifying hypothesis was that $x$ would be squarefree. He referred to these problems with the name Ozanam, from somewhere in Dickson's History.
The reason something can be done is that it is possible to regard odd/even exponents as numbers $0,1$ and perform row reductions using Gauss methods in the field with two elements.
So, take the primes up to $100,$ say. Each prime gets a column in a certain matrix. For each prime, factor (in my case ) $\sigma(p^3),$ and assign a row for every prime that comes up in at least one of these factorizations.
The column for prime $p$ gets mostly $0$ entries, but gets a $1$ at every prime factor of  $\sigma(p^3)$ that occurs to an odd power. 
An answer to the original problem is a vector in the nullspace of this matrix, that is a vector consisting of either $0$ or $1$ that is sent to all $0$ over the field with two elements. This is, simply, a choice of the prime factors of a (squarefree) solution.
None of this was easy, but the program was lightning fast. 
One could customize this for $\sigma(x)$ = y^2 by dropping the $x^3$ aspect, instead, for a few primes of interest, replace the column for $p$ by the column for the preferred $p^a.$ I did some of that at the time. It is less elegant, there are some nice features to keeping homogeneity of exponents. 
A: Here is where I approached a more general problem from an algorithm standpoint. However, the discussion involved should answer your question.
