I am trying to write a python program to solve for all the positive integer solutions of a polynomial of the form x_1 + 2x_2 + ... + Tx_T = a subject to x_0 + x_1 + ... + x_T = b.
For example, by brute-force I can show that the only positive integer solution to x_1 + 2x_2 = 4 subject to x_0 + x_1 + x_2 = 2 is (0,0,2). However, I would like a program to do this on its own.
Any suggestions?
You can use Numpy Linear Algebra to solve a system of equations, the least-squares solution to a linear matrix equation. In your case, you can consider the following vectors:
import numpy as np
# range(T): coefficients of the first equation
# np.ones(T): only 'ones' as the coefficients of the second equation
A = np.array([range(T), np.ones(T)) # Coefficient matrix
B = np.array([a, b]) # Ordinate or “dependent variable” values
And then find the solution as follows:
x = np.linalg.lstsq(A, B)[0]
Finally you can implement a solve method passing the variables T, a b:
import numpy as np
def solve(T, a, b):
A = np.array([range(T), np.ones(T)])
B = np.array([a, b])
return np.linalg.lstsq(A, B)[0]
Integer solutions: If you want only integer solutions then you are looking at a system of linear diophantine equations.
every system of linear Diophantine equations may be written: AX = C, where A is an m×n matrix of integers, X is an n×1 column matrix of unknowns and C is an m×1 column matrix of integers.
Every system of linear diophantine equation may be solved by computing the Smith normal form of its matrix.
(0, 0, 2)right?