What other language allows you to do Rational matrix inversion?

Specific Domains Numerics
1

I was just doing some Hackerrank for fun, in particular this one.

Basically, you can describe the problem as a sequence of recurrence relations, and they told you that the solution is in the form of a rational number. So far so good.

TLDR: If M and J are Matrix{Rational{Int}} then M\J is also Matrix{Rational{Int}}! I found this amazing! Is Julia fairly unique in this aspect? I know Mathematica can probably do this too but I wonder if there are any others? This makes Julia really powerful

Long version: One way to solve a recurrence relation is of course to use matrices. So I build a Matrix{Rational{Int}} that represents the recurrence relation and I just inverted it. The awesome thing was the inverse of Matrix{Rational{Int}} is also Matrix{Rational{Int}} and I just solved the whole system by doing M\J where M is the recurrence relations and J is the RHS of the recurrence relations and the solution of M\J is also Matrix{Rational{Int}}! Amazing!

Unfortunately, Hackerrank doesn’t accept Julia solutions for this problem, but I’ve posted my solution in the comment.

2 Likes
2

Apparently they support Julia for some questions: have you asked them about adding Julia support for these ones?

3

I have emailed them. And I can see Julia for many questions (running v0.6 though), but just not for this.

Also they don’t allow the use of packages.

4

They updated the code with Julia support! And v1 as well. Now I need check why my solution doesn’t work!

1 Like
5

Sagemath is really good and relatively performant for this kind of thing. It has all of the python drawbacks, though, and is afaik still stuck on 2.7 (good and fast external libraries, but two-language problem for anything you want to be fast).

1 Like
6

Mathematica does it.

1 Like
7

@Juan yeah I already mentioned Mathematica in my post. I just don’t know the syntax. It would be great if you can post some code to do M\J in Mathematica where M and J are Matrices with Rational number and of course M\backslash J = M^{-1}J.

Looks like Matlab has symbolic maths so it’s doable as well. Just not sure how easy the syntax is?

8

I’m not an expert but just:

mat = {{1, 2, 3}, {4, 5, 6}, {7, 8, 10}}

Inverse[mat]

{{-(2/3), -(4/3), 1}, {-(2/3), 11/3, -2}, {1, -2, 1}}

add // MatrixForm on the right if you want a formatted matrix.

Is this what you mean?
Anyway you’d better asking Mathematica questions on a dedicated forum.

1 Like
9

Here’s mathematica:

m = {{ 5/4,  4/3,  4/5},
   {4/5,  5/1,  1/1},
   {1/5,  7/1,  3/4}}
j = {{ 10/3,   1/1,  8/3},
   { 5/1,   5/9,  1/5},
   { 1/10 , 2/3,  1/1 }}
LinearSolve[m, j]

which returns this:

{{-(2040/157), 8900/3297, 12496/1099}, {-(2137/785), 1349/3297, 9127/
  5495}, {4554/157, -(36100/9891), -(18904/1099)}}

Julia, from which I got those numbers!

julia> M = rand(1:10, 3,3) .// rand(1:10, 3,3)
3×3 Array{Rational{Int64},2}:
 5//4  4//3  4//5
 4//5  5//1  1//1
 1//5  7//1  3//4

julia> J =rand(1:10, 3,3) .// rand(1:10, 3,3)
3×3 Array{Rational{Int64},2}:
 10//3   1//1  8//3
  5//1   5//9  1//5
  1//10  2//3  1//1

julia> M \ J
3×3 Array{Rational{Int64},2}:
 -2040//157    8900//3297   12496//1099
 -2137//785    1349//3297    9127//5495
  4554//157  -36100//9891  -18904//1099

It’s very likely you could do this in Maple too, but I don’t have that right now. Here’s a random link.

10

Please be aware of easy undetected integer overflows, when you do M \ J, especially if the matrix size is more than 10 or so… In this case Matrix{BigInt} should do.

1 Like
11

I’d recommend against Matrix{Rational{BigInt}}. Generic code that falls back on BigInt is just slow in julia. For standard cases (like blas), you really want specialized libraries like the ones used by nemo.

12

Could maybe use Rational{SaferIntegers.SafeInt}, and if that errors fall back to something slower. Julia’s composability shines here I think (although it looks like @JeffreySarnoff put conscious effort into improving Rational support in SaferIntegers recently).

1 Like
13

I did end up using BigInt for my solutions

14

Another language you can use is REDUCE, for which I made a Julia interface Reduce.jl

You can try it out from within the Julia language, it supports inverting a rational matrix and also symbolic matrix, with arbitrary precision also.

julia> using Reduce
Reduce (Free CSL version, revision 4590), 11-May-18 ...

julia> [:x :y; :z :w] |> Reduce.Algebra.inv
2×2 Array{Any,2}:
 :(w / (w * x - y * z))   :(-y / (w * x - y * z))
 :(-z / (w * x - y * z))  :(x / (w * x - y * z))

The Reduce.Algebra.inv method should default to Julia’s internal methods for non-symbolic input, but if you really want to ask REDUCE to invert your rational matrix for you from Julia, you can do it by converting the Matrix to RExpr first

julia> M = rand(1:10, 3,3) .// rand(1:10, 3,3);

julia> RExpr(M)

[ 1          1 ]
[---   1    ---]
[ 3          5 ]
[              ]
[ 7    10      ]
[---  ----   2 ]
[ 3    9       ]
[              ]
[      7       ]
[ 1   ----   9 ]
[      10      ]


julia> string(ans)
"mat((1/3,1/1,1/5),(7/3,10/9,2/1),(1/1,7/10,9/1))"

which is internally represented as REDUCE source code for a rational matrix

Then, you need to enable the Reduce.Rational(true) switch, so that the Reduce parser knows you want to parse Rational things instead of Floats for conversions

julia> Reduce.Rational(true);

julia> RExpr(M) |> Reduce.Algebra.inv

[  - 3870     3987       - 800 ]
[---------   ------    --------]
[  7213       7213       7213  ]
[                              ]
[  8550       - 1260      90   ]
[ ------    ---------   ------ ]
[  7213       7213       7213  ]
[                              ]
[  - 235      - 345     2650   ]
[--------   --------   ------- ]
[  7213       7213      21639  ]

julia> eval(parse(ans))
3×3 Array{Rational{Int64},2}:
 -3870//7213   3987//7213  -800//7213 
  8550//7213  -1260//7213    90//7213 
  -235//7213   -345//7213  2650//21639

There you go, that is how you would do it with the REDUCE language, from within Julia.

It supports arbitrary precision, but due to the parser interface is limited in performance.

15

One of my lecturer taught Reduce for a course. I forgot all about it!

16

tis true

note: ops over a matrix of Rational{SafeInt} values will check for overflow (which can happen as Rationals tend to grow in digits).

I suggest using Rational{SafeInt128}s to hold values that are Rational{SafeInt16}s, and operate on those to sess out algorithmic design hitches.

q16 = Rational{SafeInt16}(3, 5) 
#       == SafeInt16(3)//SafeInt16(5)
q128 = Rational{SafeInt128}(q16)

Broadcasting the size conversion works, too.