I am personally content with the following solution: I calledsage in Python and letsage compute what I want.sage, being math-oriented, is rather versatile in computations involving fields other than reals.

Below is my scriptcompute_intarrs.py and it requiressage be installed. Be aware it is a little slow.

import subprocessimport reimport numpy as np# construct a numpy matrixmat = np.matrix([[1,-1],[-1,1]])# convert the matrix into a string recognizable by sagematstr = re.sub('\s|[a-z]|\(|\)', '', mat.__repr__())# write a (sage) python script "mat.py";# for more info of the sage commands: # www.sagemath.org/doc/faq/faq-usage.html#how-do-i-import-sage-into-a-python-script# www.sagemath.org/doc/tutorial/tour_linalg.htmlf = open('mat.py', 'w')f.write('from sage.all import *\n\n')f.write('A = matrix(ZZ, %s)\n\n' % matstr)f.write('print A.kernel()')  # this returns the left nullspace vectorsf.close()# call sage and run mat.pyp = subprocess.Popen(['sage', '-python', 'mat.py'], stdout=subprocess.PIPE)# process the output from sagearrstrs = p.communicate()[0].split('\n')[2:-1]arrs = [np.array(eval(re.sub('(?<=\d)\s*(?=\d|-)', ',', arrstr)))         for arrstr in arrstrs]print arrs

Result:

In [1]: %run compute_intarrs.py

[array([1, 1])]

You can do some pretty cool things withdtype = object and thefractions.Fraction class, e.g.

>>> A = np.array([fractions.Fraction(1, j) for j in xrange(1, 13)]).reshape(3, 4)>>> Aarray([[1, 1/2, 1/3, 1/4],       [1/5, 1/6, 1/7, 1/8],       [1/9, 1/10, 1/11, 1/12]], dtype=object)>>> B = np.array([fractions.Fraction(1, j) for j in xrange(1, 13)]).reshape(4, 3)>>> Barray([[1, 1/2, 1/3],       [1/4, 1/5, 1/6],       [1/7, 1/8, 1/9],       [1/10, 1/11, 1/12]], dtype=object)>>> np.dot(A, B)array([[503/420, 877/1320, 205/432],       [3229/11760, 751/4620, 1217/10080],       [1091/6930, 1871/19800, 1681/23760]], dtype=object)

Unfortunately thenp.linalg module converts everything tofloat before doing anything, so you can't expect to get solutions directly as integers or rationals. But you can always do the following after your computations:

def scale_to_int(x) :    fracs = [fractions.Fraction(j) for j in x.ravel()]    denominators = [j.denominator for j in fracs]    lcm = reduce(lambda a, b: max(a, b) / fractions.gcd(a, b) * min(a, b),                 denominators)    fracs = map(lambda x : lcm * x, fracs)    gcd = reduce(lambda a, b: fractions.gcd(a, b), fracs)    fracs = map(lambda x: x / gcd, fracs)    return np.array(fracs).reshape(x.shape)

It will be slow, and very sensitive to round-off errors:

>>> scale_to_int(np.linspace(0, 1, 5)) # [0, 0.25, 0.5, 0.75, 1]array([0, 1, 2, 3, 4], dtype=object)>>> scale_to_int(np.linspace(0, 1, 4)) # [0, 0.33333333, 0.66666667, 1]array([0, 6004799503160661, 12009599006321322, 18014398509481984], dtype=object)

You could mitigate some of that using thelimit_denominator method ofFraction, but probably will not be all that robust.

• python
• numpy
• integer
• scipy
• eigenvector