I have a particular system of linear Diophantine equations in $n$ variables for which I need to find all nonnegative integer solutions.
Specifically, they are Frobenius equations, meaning the coefficients are natural numbers and nonnegative integer solutions. The equations are:
$$\left(\begin{matrix} 1 & 1 & ... & 1\\ 1 & 2 & ... & n \\ \end{matrix}\right) \left(\begin{matrix} x_1 \\ \vdots \\ x_n \end{matrix}\right) =\left(\begin{matrix} M-r\\ M \end{matrix}\right) $$
n and M are fixed, and I'll need to solve for each of r=0,..., floor$\left(\frac{M(n-1)}{n}\right)$.
Ideally, I will be able to solve relatively large $M$ (as large as computation and storage allows).
I had some luck using Mathematica's built-in Frobenius equation solver to generate all solutions for the second equation and then pulling out solutions satisfying $x_1 + ... + x_n=M-r$ for each $r$. For $n=10$ (which is sufficient for my needs), I was able to solve for as high as $M=100$, but going much larger made computation time (and space) troublesome.
Is there a way to take advantage of the existence of the second equation, linearity, or some other detail that isn't being considered to solve this more efficiently?
Is solving it with $n=10, M=1000$ a pipe dream on a strong desktop computer? Would it even be attainable on more substantial hardware?
Thanks!