Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Welcome To Ask or Share your Answers For Others

Categories

0 votes
420 views
in Technique[技术] by (71.8m points)

numpy - np.linalg.det() different result than my calculation

I am trying to understand how computing the determinant works. I have:

b = np.array([ [1,2,3], [2,3,4], [5 6,7] ])
np.linalg.det(b)
>>> 1.7763568394002396e-15

Now, this is how I understood it:

>>>>>>> array([[1, 2, 3],
               [2, 3, 4],
               [5, 6, 7]])
>>>>>>> |b| = |1  2  3|   = 1|3  4| - 2|2  4| + 3|2  3|
              |2  3  4|      |6  7|    |5  7|    |5  6|
              |5  6  7|
        |b|               = 1*(3*7 - 4*6) - 2*(2*7 - 4*5) + 3*(2*6 - 3*5)
        |b|               = 1*(   - 3   ) - 2*(    - 6  ) + 3*(    - 3  )    
        |b|               =   (   - 3   ) -   (    - 12 )   + (    - 9  )
        |b|               =       -3               +12             - 9     = 0      

Haven't I understood it correctly?

Why det() is returning 1.77 while I am ending up with 0 ?

Should it be maybe the following instead?

>>>>>>> |b| = |1  2  5|   = 1|3  6| - 2|2  6| + 5|2  3|
              |2  3  6|      |4  7|    |3  7|    |3  4|
              |3  4  7|
question from:https://stackoverflow.com/questions/65937458/np-linalg-det-different-result-than-my-calculation

与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
Welcome To Ask or Share your Answers For Others

1 Reply

0 votes
by (71.8m points)

What the above comments meant is that numpy output:

b = np.array([ [1,2,3], [2,3,4], [5,6,7] ], dtype=np.float64)
np.linalg.det(b)
>>> 1.7763568394002396e-15

may roughly be seen as 0. +/- roundoff_errors where roundoff_errors are pretty close to np.finfo(np.float64).eps = 2.22e-16 which corresponds to machine accuracy for float64 datatype. Indeed, np.linalg.det does not compute determinant analytically but relies on a numerical approximation using LU factorization of b (as stated in the documentation). I don't have looked at the source file to see what numpy exactly does but it surely boils down to something like:

from scipy.linalg import lu_factor
b_LU = lu_factor(b)

def lu_det(A_LU):

    # The determinant of a permutation matrix is (-1)**n where n is the number of permutations.
    # Recall a permutation matrix is a matrix with a one in each row, and column, and zeros everywhere else.
    # => if the elements of the diagonal are not 1, they will be zero, and then there has been a swap
    nswaps = len(A_LU[1]) - sum(A_LU[1]==np.arange(len(A_LU[1])))
    detP = (-1)**nswaps

    # The determinant of a triangular matrix is the product of the elements on the diagonal.
    detL = 1 # The unit diagonal elements of L are not stored.
    detU = np.prod(np.diag(A_LU[0]))
    # Matrix containing U in its upper triangle, and L in its lower triangle.

    # The determinant of a product of matrices is equal to the product of the determinant of the matrices.
    return detP * detL * detU

lu_det(b_LU)
>>> 4.4408920985006196e-17

which again is non-zero but quite close indeed! If you want a true zero, you may want to have a look at sympy for symbolic calculations, which gives simply:

import sympy as sp

sp.Matrix(b).det()
>>> 0

Hope it is clearer to you now.


与恶龙缠斗过久,自身亦成为恶龙;凝视深渊过久,深渊将回以凝视…
OGeek|极客中国-欢迎来到极客的世界,一个免费开放的程序员编程交流平台!开放,进步,分享!让技术改变生活,让极客改变未来! Welcome to OGeek Q&A Community for programmer and developer-Open, Learning and Share
Click Here to Ask a Question

...