CouscousPie
CouscousPie

Reputation: 5

Matrix multiplied by its Inverse doesn't return Identymatrix (but the other way around does)

So, I'm just figuring things out and can't wrap my head around the fact that the following code doesn't give me back an Identymatrix.

import numpy as np

b = np.array([[2,3,5],[1,2,1],[4,6,3]])
bm = np.asmatrix(b)
print(bm)

c = np.linalg.inv(bm)
cm = np.asmatrix(c)
print(cm)

bm*cm

What baffles me is that

cm*bm

gives back the identymatrix. But why? For matrix multiplication the opposite had to be true

A*A^(-1) = I

Where A is bm and A^(-1) is cm.

Upvotes: 0

Views: 913

Answers (2)

hpaulj
hpaulj

Reputation: 231385

In [756]: b = np.array([[2,3,5],[1,2,1],[4,6,3]])
In [757]: b
Out[757]: 
array([[2, 3, 5],
       [1, 2, 1],
       [4, 6, 3]])

Sticking with the np.array (np.matrix use is discouraged these days)

In [758]: c = np.linalg.inv(b)
In [759]: c
Out[759]: 
array([[ 0.        , -3.        ,  1.        ],
       [-0.14285714,  2.        , -0.42857143],
       [ 0.28571429,  0.        , -0.14285714]])

Note that the values are floats, which even with float64 have limited precision.

Matrix multiplication of the two. (@ is matrix multiplication for both ndarray and np.matrix)

In [760]: b@c
Out[760]: 
array([[ 1.00000000e+00,  0.00000000e+00,  0.00000000e+00],
       [ 0.00000000e+00,  1.00000000e+00, -5.55111512e-17],
       [ 0.00000000e+00,  0.00000000e+00,  1.00000000e+00]])
In [761]: np.allclose(b@c,np.eye(3))
Out[761]: True

allclose is recommended when testing equality with floats. It allows for differences like -5e-17.

In sympy which uses symbols, algebra, and rationals (fractions):

In [38]: B = Matrix([[2,3,5],[1,2,1],[4,6,3]])

In [39]: B
Out[39]: 
⎡2  3  5⎤
⎢       ⎥
⎢1  2  1⎥
⎢       ⎥
⎣4  6  3⎦

In [40]: B.inv()
Out[40]: 
⎡ 0    -3   1  ⎤
⎢              ⎥
⎢-1/7  2   -3/7⎥
⎢              ⎥
⎣2/7   0   -1/7⎦

In [41]: B * B.inv()
Out[41]: 
⎡1  0  0⎤
⎢       ⎥
⎢0  1  0⎥
⎢       ⎥
⎣0  0  1⎦

Note those 1/7 fractions. Scaling the float inverse:

In [762]: c*7
Out[762]: 
array([[  0., -21.,   7.],
       [ -1.,  14.,  -3.],
       [  2.,   0.,  -1.]])

And doing integer matrix multiplication:

In [765]: b@((c*7).astype(int))
Out[765]: 
array([[7, 0, 0],
       [0, 7, 0],
       [0, 0, 7]])

We can scale that back to 1:

In [767]: b@((c*7).astype(int))/7
Out[767]: 
array([[1., 0., 0.],
       [0., 1., 0.],
       [0., 0., 1.]])

But if we scale before the multiplication:

In [768]: b@((c*7).astype(int)/7)
Out[768]: 
array([[1.00000000e+00, 0.00000000e+00, 2.22044605e-16],
       [0.00000000e+00, 1.00000000e+00, 5.55111512e-17],
       [0.00000000e+00, 0.00000000e+00, 1.00000000e+00]])

Upvotes: 0

roxolid
roxolid

Reputation: 93

A*A^(-1) = I and A^(-1)*A = I

should both be true.

I get somethin like this for the first multiplication:

[ 1.00000000e+00,  0.00000000e+00, -5.55111512e-17],
        [ 0.00000000e+00,  1.00000000e+00, -5.55111512e-17],
        [ 0.00000000e+00,  0.00000000e+00,  1.00000000e+00]])

Seems to be a numeric problem, but -5*10^-17 is preatty close to 0, so if you were to round your results, you would get I.

By the way, you dont need to convert c to a Matrix as it allready is.

Upvotes: 2

Related Questions