Matrix inversion modulo p; Linear systems modulo p
Vector 18:1, page 134
Transcribed from the printed issue; not yet reviewed. Read it in the PDF of the issue, page 134
Summary¶
The two problems announced in the title are linked. To solve the one is to solve the other. On the other hand solving the first leads to the resolution of another problem that had not yet been solved [1] before this present attempt. This article was first published in Les Nouvelles d’APL N° 35.
Let A be a matrix of integers of order n, p an integer and V a vector of integers of n elements all less than p. Given W = p|A+.###V, is there a matrix B such that p|B+.###W restores the original vector V?
The answer is yes, if A is well formed (but not in the way one normally thinks – see [2]). We give here a recipe for the calculation of B under this assumption. We will call B the inverse modulo p of A (the term inverse of A modulo p would be more ambiguous).
Introduction¶
Such a system can be solved manually up to order 3. For any order n, despite reading through a large bibliography and searching the web, nowhere have we found a solution to this problem. As usual, using APL and J tools (that, a long time ago, I christened ‘faitout-vite’ in French, possibly SpeedyGee in English), I found the general solution in one day, first of all in APL. As the numbers grow very quickly, eventually I had to resort to J, which lets one work in extended precision, that is integers without limit. As they say, if you want it done, do it yourself, and I think we should not hesitate when we have tools such as APL and J at our disposal.
The work here presented served as preparation for a paper where the techniques developed are put to a practical application in the exciting field of cryptology. [2]
We will deal in turn first with the APL solution and then the J solution.
The APL solution leads to a deadlock when the matrix A to invert (modulo p) is of a order high enough, because we would need numbers larger than the limit of APL integers. However, we will not pass in silence over this work we had to do.
1. APL solution¶
As we need the integer matrix determinant and inversion, we pick up our dusty old notebooks and write these functions in APL.
Calculation of the adjoint matrix:
Calculation of the cofactors:
∇ CO←COFACTORS A;I;J;⎕IO
[1] CO←⍳0 ⋄ I←1
[2] INC:J←0
[3] IT:J←J+1
[4] CO←CO,A COFACTOR I,J
[5] →(J<(⍴A)[1])/IT
[6] I←I+1
[7] →(I≤1↑⍴A)/INC
[8] CO←(⍴A)⍴CO
[9]
∇
Calculation of a cofactor:
Calculation of the determinant:
Calculation of a minor:
Calculation of the inverse:
Thus we are well prepared and are sure that we a re working with integers.
We can now write the function that calculates the inverse modulo p of a matrix.
Calculation of the inverse matrix modulo p¶
Taking into account successive vectors from the first to the last line of the unit matrix of order n, we have the next solution. We call it DIO in honour of the Greek mathematician Diophantus.
∇ Z←M DIO P;⎕IO;I;T;U;N;IT;MU;D
[1] ⎕IO←1⋄ T←⍉⌽M ⋄N←1↑⍴M
[2] U←(⍳N)∘.=⍳N ⋄ IT←INV T ⋄D←DET M
[3] MU←⊖U ⋄ Z←(N,N)⍴0 ⋄ I←1
[4] INC:Z[I;]←(P|D×⍳P)⍳⌊(P|((IT)+.×MU[I;])×D)
[5] →(N≥I←I+1)/INC
[6]
∇
(This is the function I searched for everywhere and eventually I had to write for myself as well.)
Test of the function¶
Let us choose p=256.
We take a random integer matrix of order 3, whose elements are in the range -255 to 255
Let us calculate the matrix B inverse modulo 256 of A
Let us take a random vector V whose elements are less than 256.
Let us apply to it P|A+.×V
The goal of the game is to recover V by applying P|B+.×W.
Success!
Do not rejoice too soon!
Unfortunately already at order 6, the calculations will overflow, causing the inverse to be wrong.
⎕←A←6 6⍴36?256
33 192 116 135 55 11
173 173 238 97 132 211
8 13 135 171 1 97
17 106 175 150 237 215
134 23 166 106 178 232
194 66 12 187 83 161
A DIO P
60 156 224 257 136 257
11 257 202 124 257 257
257 109 257 55 257 257
54 71 257 42 68 173
229 141 77 202 257 67
148 257 257 14 77 257
Let us start again with another (because a single element of 257 is enough to make the calculation of b wrong).
A←6 6⍴36?255
131 214 7 108 184 118
76 83 224 112 96 82
124 105 100 135 204 119
64 230 248 113 159 225
176 207 142 103 92 89
15 52 110 115 78 67
A DIO P
257 257 257 257 257 257
257 257 257 257 257 257
257 257 257 257 257 257
257 257 110 257 29 257
257 257 257 257 257 257
257 257 79 257 45 257
It’s about time to gear up and switch to J, the only one which can save us, since it can store and perform calculations on integers as large as we need.
2. J Solution¶
Let us rewrite our routines.
For the determinant, no sweat:
Luckily for us, this train is special-cased in the latest releases of the interpreter produced by Jsoftware and we don’t have to hide the beauty of the mathematical definition behind a thousand lines of efficient, but ugly, code.
J does not have among its primitives the equivalent of APL’s vertical mirror, so let us re-write it in tacit form as:
We will trust the inverse matrix calculated by the J primitive %.
dio=: 3 : 0
:
t =: |: mv x.
n =:0{$x.
u =: (i.n) =/ i.n
mu=:|.u
it=: %. t
d =: det x.
b =:(n,n)$0
i =:0
while. i<n do.
b =: (y.|d*i.y.) i. (y. |((it)+/ .*(i{mu))*d) (i}) b
i =:i+1
continue.
end.
)
Now we can try our teeth with the first matrix that did not work in APL.
a
33 192 116 135 55 11
173 173 238 97 132 211
8 13 135 171 1 97
17 106 175 150 237 215
134 23 166 106 178 232
194 66 12 187 83 161
a:=x: a NB. Conversion into integers without limit.
p =:256x
b =:a dio p
b
175 50 54 130 174 65
46 166 248 202 13 54
198 72 47 82 9 185
121 146 33 53 87 179
123 139 24 136 171 182
186 153 237 33 96 228
Let us pick a random vector v:
Let us calculate w and verify that we get v back:
In order to show that it works on matrices of high order, just to whet our appetite for the paper on cryptology, let us take a matrix of order 10 and do a little test.
+a=:10 10$100?256
37 149 215 150 243 141 37 250 104 36
144 64 124 118 245 32 50 81 160 32
166 158 204 63 121 99 51 7 229 108
36 241 104 33 225 23 41 18 93 64
34 199 116 89 115 206 237 166 54 173
231 63 219 120 129 153 208 192 117 242
161 112 210 175 179 251 243 217 73 137
131 26 105 147 223 112 186 221 182 204
180 189 4 225 133 118 16 181 124 170
173 50 233 220 226 138 35 114 252 54
+b=:a dio p
168 146 70 169 227 106 107 63 184 81
38 129 50 213 66 33 80 54 219 157
42 75 118 174 151 243 131 72 126 140
44 206 38 39 196 223 76 214 6 241
74 98 77 108 76 137 198 112 215 147
87 209 5 217 234 46 96 214 236 223
80 31 80 237 212 111 6 81 40 234
78 177 51 40 80 39 202 42 145 55
35 123 194 69 112 0 144 149 2 18
60 12 17 23 13 222 78 51 166 227
+v=:10?256
234 79 158 41 220 6 242 228 125 56
+w=:p|a+/ . *v
241 162 174 217 22 28 170 135 124 136
p|b+/ . *w
234 79 158 41 220 6 242 228 125 56
NB. In both matrices, one can replace any number n (especially those that are greater than 128) by the negative number n - 256.
Let us do it on a line of a taken randomly.
u =: _25 63 _37 120 _127 _103 _48 _64 117 _14
u =: x: u
a =: u (5}) a
a
37 149 215 150 243 141 37 250 104 36
144 64 124 118 245 32 50 81 160 32
166 158 204 63 121 99 51 7 229 108
36 241 104 33 225 23 41 18 93 64
34 199 116 89 115 206 237 166 54 173
_25 63 _37 120 _127 _103 _48 _64 117 _14
161 112 210 175 179 251 243 217 73 137
131 26 105 147 223 112 186 221 182 204
180 189 4 225 133 118 16 181 124 170
173 50 233 220 226 138 35 114 252 54
1) The transformed matrix b is identical:
a dio p
168 146 70 169 227 106 107 63 184 81
38 129 50 213 66 33 80 54 219 157
42 75 118 174 151 243 131 72 126 140
44 206 38 39 196 223 76 214 6 241
74 98 77 108 76 137 198 112 215 147
87 209 5 217 234 46 96 214 236 223
80 31 80 237 212 111 6 81 40 234
78 177 51 40 80 39 202 42 145 55
35 123 194 69 112 0 144 149 2 18
60 12 17 23 13 222 78 51 166 227
2) the transformed vector w is the same:
We shall not pass over in silence the ulterior motive that some may have had reading this paper…
It’s no accident that we have chosen p=256 because the residue modulo 256 is a number in the range 0 to 255 and can represent the index in origin 0 of any character of a computer in the ⎕AV.
If b is the inverse modulo p of a, it is necessary to know that it is symmetrical: then a is the inverse modulo p (the same p of course) of b.
Let us verify this by working on the same vector v again but let us use b instead of a as the transformation matrix. Except in case of an unlikely coincidence, we should not have the same w.
The w transformed with a was:
It does not matter if the w calculated is not the same, we are going to find our v by the magical operation:
One concludes immediately that the matrices a and b can be kept secret separately. Person B for example, who knows how to invert modulo p, sends to his correspondent (A) the matrix a that he can use to transform messages for B’s eyes only. B reconstitutes the messages with b but can send messages to A by using b. A will then reconstitute the message with a and does not need to know b.
Of course what is sent down the wire will not be the numeric vector w but the characters situated at the indices contained in w in the ⎕AV (or maybe something cleverer – see [2]). In other words, the person A would transmit the message T[w] (reads naturally T indexed by w), where T is the atomic vector in our example.
A last paragraph on an interesting property of the transformation: of course one could take matrices a whose elements are greater than -255 and less than 255 but it would be as if he took the matrix 256|a.
References¶
[1] Michel J. Dumontier. Inversion des matrices Booléennes. Les Nouvelles d’APL N° 35 p.55. Translated as Inversion of Boolean Matrices, Vector 18.1
[2] Michel J. Dumontier. Un nouveau chiffrement court échappant à la cryptanalyse. Les Nouvelles d’APL N° 35 p.65. Translated as MUCIPHER: A new short cipher that escapes cryptanalysis, Vector 18.2.