I'm trying find the minimum sum of $x_{1} + x_{2} + ... + x_{n}$ where these are a solution to a linear system of two equations.
System of linear equations in general form: $$ a_{11}x_{1} + a_{12}x_{2} + ... a1_{n}x_{n} = b_{1} \\ a_{21}x_{1} + a_{22}x_{2} + ... a2_{n}x_{n} = b_{2} $$ The unknown variables $(x_{i})$ and coefficients $(a_{ji})$ are non-negative whole numbers, i.e:
$$ x_{i}\in\mathbb{Z}_{\geq0} \\ a_{ij}\in\mathbb{Z}_{\geq0} $$ As for the constant terms $(b_{j})$, they are positive whole numbers, i.e: $$ b_{j}\in\mathbb{Z}_{\gt0} $$
The number of unknown variables may be any positive number $(n\in\mathbb{Z}_{\gt0})$.
Since each unknown variable is non-negative, there are a finite number of solutions, which may be zero in the case the system is not solvable.
Example: $$ 10x_{1} + 5x_{2} + 0x_{3} = 10 \\ 2x_{1} + 0x_{2} + 1x_{3} = 2 $$ The only two solutions are $$ 1: x_{1} = 1, x_{2} = 0, x_{3} = 0 \\ 2: x_{1} = 0, x_{2} = 2, x_{3} = 2 $$ where the sum of $x_{1} + x_{2} + x_{3}$ in the first solution is $1$, and $2$ in the second, making the first solution the answer.
How do I find the minimum sum of $x_{1} + x_{2} + ... + x_{n}$ that satisfy two equations with the constraints above? Solution sets are a challenge to describe when using so many unknown variables, is there a better way to solve my question? Can Matrix multiplication with Moore-Penrose pseudo-inverse be used? Any help is much appreciated!
Your problem is what's called an integer linear program (ILP), you can read about it online. In general, ILP is very hard to solve (NP-hard, to be exact) although certain cases have fast solutions. Unfortunately yours is not one of those easy cases. Just to find a solution $x$ for one of the equations, (let alone the one with minimal sum), you have to find an integer combination of the $a_{ij}$ that add up to $b_i$. This is a variation of the "subset sum" problem and also "knapsack problem," known to be NP-complete.
If you read about the knapsack problem, you will see that you can use "dynamic programming" to solve it reasonably fast if the $a_{ij}$ and $b_j$ are not too large. In your case you have two equations instead of one, but you can still take a dynamic programming approach.
NOTE: You can also use a standard ILP solver, which may work faster because they usually have a lot of heurisitics (although still guaranteed to find an optimal solution). Cplex is a commercial LP and ILP solver (although I think the free versions have limits on how many variales you can have). If nothing else, you can code your own branch-and-bound ILP solver if you don't find one online to suit your needs. (Branch-and-bound ILP solving is easy to read about online.)