-1

From some other information like: How to handle big integers in C I know that GMP is a package that enables me to handle big integers.

However, I am considering 500×500 matrix, each element of which is of length 1000 bits. Can someone let me know which package in C or Python can allow me to compute the matrix inverse?

Community
  • 1
  • 1
user4478
  • 1,289
  • 2
  • 11
  • 13
  • I would use MatLab for this. – Jordan Kaye Jan 11 '13 at 18:48
  • I also used matlab to verify the answer, but I currently need some C or Python implementation – user4478 Jan 11 '13 at 19:06
  • 1
    This question needs a better definition of your problem and intended solution. Answers could be as simple as "compute a IEEE 754 floating point approx of your matrix and use a standard linear algebra package" up to some esoteric approach: it all depends on the property of your matrix and the errors in the computed solution you are willing to accept. – Stefano M Jan 14 '13 at 14:26

2 Answers2

1

As casevh pointed out, this is highly dependent on whether you are want an approximate answer or an exact answer.

C libraries that might be suited for computing an exact inverse of a matrix of this size include FLINT, PARI and IML. FLINT takes this long to invert a random n by n integer matrix with 1000-bit entries:

n = 13: 1.6 seconds

n = 62: 37 seconds

n = 125: 858 seconds

So n = 500 should take about 120 hours.

Fredrik Johansson
  • 25,490
  • 3
  • 25
  • 17
0

I'm not sure from your question if you are trying to calculate the inverse using floating point numbers with 1000 bits of precision, or if you want to use fractions.

For the former, I would try the combination of mpmath and gmpy. The following example creates a random 3x3 matrix with 100 bits of precision and calculates the inverse.

>>> from mpmath import *
>>> mp.prec = 100
>>> a = randmatrix(3)
>>> a**-1
matrix(
[['-2.9551532644739824344170194538', '-2.30592481828272627130234272', '4.618043964228917637573009707'],
 ['12.025269724235234315848646394', '3.3570351961866008157001332066', '-10.59068013761452569858474265'],
 ['-6.672524995207041946321450867', '-0.57061969261482431665347164675', '5.508560423258977568243759022']])
>>> a * a**-1
matrix(
[['1.0', '-3.8920288101260775500720697474e-34', '5.467369412738327810686299975e-34'],
 ['-2.0267823641769153744296258652e-33', '1.0', '-1.4939276515262026082238333732e-33'],
 ['-2.1979577599335336964264104571e-33', '-6.0264087812469052430270713117e-34', '1.0']])
>>>

Note: mpmath 0.17 will only work with gmpy. To use gmpy2 (the next major release), you will need to use the source repository of mpmath.

If you want an exact rational inverse, you can try sympy. I would be concerned that the size of the fractions will increase and it will be much slower than the floating point inverse. sympy uses mpmath and gmpy behind the scenes.

casevh
  • 11,093
  • 1
  • 24
  • 35