Solution 1:

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.