Solve system of linear integer equations in Python

I am looking for a method to solve a system of linear equations in Python. In particular, I am looking for the smallest integer vector that is larger than all zeros and solves the given equation. For example, I have the following equation:

and want to solve .

In this case, the smallest integer vector that solves this equation is .

However, how can I determine this solution automatically? If I use scipy.optimize.nnls, like

A = np.array([[1,-1,0],[0,2,-1],[2,0,-1]])
b = np.array([0,0,0])
nnls(A,b)

the result is (array([ 0., 0., 0.]), 0.0). Which is also correct, but not the desired solution...

Edit: I apologize for being imprecise in certain aspects. If anyone is interested in the details, the problem comes from the paper "Static Scheduling of Synchronous Data Flow Programs for Digital Signal Processing", Edward A. Lee and David G. Messerschmitt, IEEE Transactions on Computers, Vol. C-36, No. 1, pp. 24-35, January, 1987.

Theorem 2 says

For a connected SDF graph with s nodes and topology matrix A and with rank(A)=s-2, we can find a positive integer vector b != 0 such that Ab = 0 where 0 is the zero vector.

Directly after the proof of Theorem 2 they say

It may be desirable to solve for the smallest positive integer vector in the nullspace. To do this, reduce each rational entry in u' so that its numerator and denominator are relatively prime. Euclid's algorithm will work for this.

Solution 1:

To find the exact solution that you want, numpy and scipy are probably not the best tools. Their algorithms generally work in floating point, and are not guaranteed to give the exact answer.

You can use sympy to get the exact answer to this problem. The Matrix class in sympy provides the method nullspace() that returns a list of basis vectors for the nullspace. Here's an example:

In [20]: from sympy import Matrix, lcm

In [21]: A = Matrix([[1, -1, 0], [0, 2, -1], [2, 0, -1]])

Get the vector in the nullspace. nullspace() returns a list; this code assumes that the rank of A is 2, so the list has length one:

In [22]: v = A.nullspace()[0]

In [23]: v
Out[23]: 
Matrix([
[1/2],
[1/2],
[  1]])

Find the least common multiple of the denominators of the values in v so that we can scale the result to the smallest integers:

In [24]: m = lcm([val.q for val in v])

In [25]: m
Out[25]: 2

x holds the integer answer:

In [26]: x = m*v

In [27]: x
Out[27]: 
Matrix([
[1],
[1],
[2]])

To convert that result to a numpy array of integers, you can do something like:

In [52]: np.array([int(val) for val in x])
Out[52]: array([1, 1, 2])

Solution 2:

Actually, it's just basic linear algebra.

>>> A = np.array([[1,-1,0], [0,2,-1],[2,0,-1]])    

Let's calculate the eigenvalues and eigenvector for this matrix.

>>> e = np.linalg.eig(A)
>>> e
(array([ -4.14213562e-01,  -1.05471187e-15,   2.41421356e+00]), array([[-0.26120387, -0.40824829,  0.54691816],
       [-0.36939806, -0.40824829, -0.77345908],
       [-0.89180581, -0.81649658,  0.32037724]]))    
>>>> np.round(e[0], 10)
array([-0.41421356, -0.        ,  2.41421356])

After rounding, we see that 0 is an eigenvalue of our matrix A. So, the eigenvector s for the 0-eigenvalue is a good candidate for your equation system.

>>> s = e[1][:,1]
>>> s
array([-0.40824829, -0.40824829, -0.81649658])    

Multiplying an eigenvector with the related matrix results in the eigenvector itself, scaled by the related eigenvalue. So, in the above case we see: As = 0s = 0

>>> np.round(A.dot(s), 10)
array([ 0.,  0.,  0.])

Since we are interested in an integer solution, we have to scale the solution vector.

>>> x = s / s[1]
>>> x
array([ 1.,  1.,  2.])

Hope this answer solves your problem.