Suppose $A = A^T$ and suppose the entries of $A$ are in $\mathbb{Z}^+$. I want to find all matrices $M$ with natural entries so that: $$M^2 = A$$ How can one do this? I know techniques that will get a square-root of an arbitrary matrix, but I want the full set. I want to be able to do this efficiently for large matrices ~$100 \times 100$.
Of course, the set must be finite because we are working over positive integers and the matrix $A$ gives upper bounds.
On
In a way, if you can generate one, the rest follows, since if you have one solution any other matrix $M'$ will have to adhere to the following set of equations: $$(M-M')(M+M')=M^2-M'^2 =0 $$ Now, by the definition of $M,\;M'$ they are simultaneously orthogonally diagonalizable, thus you can assume WLOG they are diagonal. even more, you know the eigen values of $M,\;M'$ sutisfy $\mu_i^2=\mu'^2_i = \lambda_i\to \mu_i=\mu'_i$ (These are the eigenvalues respectively of $M,\;M'$ and $\lambda$ is of $A$). Hence to get another matrix that answers your question is equivalent to the condition $$(\mu_i+\mu_i')\cdot(\mu_i-\mu_i')=0$$ For all of the eigenvalues. This gives a complete classification of solutions.
Find an invertible matrix $B$ and diagonal matrix $D$ such that $D=B^{-1}AB$. Then take all square roots $D_1,...,D_m$ (all of them diagonal, there are $m\le 2^n$ of these where $n$ is the size of $A$ because every non-negative number has at most 2 square roots) of $D$ and all matrices $BD_iB^{-1}$. Look which ones of these matrices are over natural numbers.