2

I have a linear system of equations like MX=N. M is a 21x21 matrix with many elements zero. When I try to solve this system with X = np.linalg.solve(M, N), it gives me this error:

numpy.linalg.linalg.LinAlgError: Singular matrix

The problem here is that the value returned by np.linalg.det(M) is 0.0. I tried two different ways to generate the M matrix and at that point I encountered a strange behavior:

i) The non-zero elements of M are calculated somewhere-else in the code. All of these elements are floats and denoted as m_1, m_2, ... , m_21. At first, I tried the following code in order to generate M:

M = np.zeros([21,21])
M[0,0] = m_1
M[0,1] = m_2
M[1,0] = m_3
M[1,4] = m_2
M[2,2] = m_2
M[2,3] = m_1
M[3,3] = m_3
M[3,5] = m_2
M[4,4] = m_4
M[4,5] = m_5
M[5,8] = m_6
M[5,13] = m_7
M[6,9] = m_6
M[6,14] = m_7
M[7,11] = m_6
M[7,12] = m_7
M[8,8] = m_8
M[8,9] = m_9
M[8,11] = m_10
M[9,6] = m_11
M[9,8] = m_12
M[9,20] = m_13
M[10,5] = m_11
M[10,10] = m_12
M[10,19] = m_13
M[11,19] = m_14
M[11,20] = m_15
M[12,8] = m_15
M[12,10] = m_14
M[13,16] = m_4
M[13,17] = m_17
M[14,7] = m_15
M[14,17] = m_16
M[15,16] = m_18
M[15,18] = m_7
M[16,17] = m_19
M[16,18] = m_20
M[17,4] = m_14
M[17,16] = m_16
M[18,11] = m_12
M[18,15] = m_13
M[19,12] = m_20
M[19,15] = m_21
M[20,7] = m_19
M[20,13] = m_20
M[20,20] = m_21

Determinant of this matrix calculated by np.linalg.det(M) is zero.

ii) Then I replaced the non-zero elements (m_1, ... , m_21) with the corresponding numeric values to see if the determinant will change. Here is the code:

 M = np.matrix([[-88.89714245, 33.72326786, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], #1
                [-139.63175129, 0, 0, 0, 33.72326786, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],#2
                [0,0,33.72326786, -88.89714245, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], #3
                [0, 0, 0, -139.63175129, 0, 33.72326786, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],#4
                [0, 0, 0, 0, 98.58344885, 55.0147276, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], #5
                [0, 0, 0, 0, 0, 0, 0, 0, 114.92510983, 0, 0, 0, 0, 66.13785145, 0, 0, 0, 0, 0, 0, 0], #6
                [0, 0, 0, 0, 0, 0, 0, 0, 0, 114.92510983, 0, 0, 0, 0, 66.13785145, 0, 0, 0, 0, 0, 0], #7
                [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 114.92510983, 66.13785145, 0, 0, 0, 0, 0, 0, 0, 0], #8
                [0, 0, 0, 0, 0, 0, 0, 0, 28.52149986, -96.35068993, 0, 67.82919006, 0, 0, 0, 0, 0, 0, 0, 0, 0], #9
                [0, 0, 0, 0, 0, 0, 83.66136319, 0, 95.15580459, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -178.81716778], #10
                [0, 0, 0, 0, 0, 83.66136319, 0, 0, 0, 0, 95.15580459, 0, 0, 0, 0, 0, 0, 0, 0, -178.81716778, 0], #11
                [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 89.26005554, 67.6481946], #12
                [0, 0, 0, 0, 0, 0, 0, 0, 67.6481946, 0, 89.26005554, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], #13
                [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,98.58344885, -153.59817645, 0, 0, 0], #14
                [0, 0, 0, 0, 0, 0, 0, 67.6481946, 0, 0, 0, 0, 0, 0, 0, 0, 0, -156.90825014, 0, 0, 0], #15
                [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -181.06296128, 0,66.13785145, 0, 0], #16
                [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -153.11049424, 35.89577791, 0, 0], #17
                [0, 0, 0, 0, 89.26005554, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -156.90825014, 0, 0, 0, 0], #18
                [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 95.15580459, 0, 0, 0, -178.81716778, 0, 0, 0, 0, 0], #19
                [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 35.89577791, 0, 0, 117.21471633, 0, 0, 0, 0, 0], #20
                [0, 0, 0, 0, 0, 0, 0, -153.11049424, 0, 0, 0, 0, 0, 35.89577791, 0, 0, 0, 0, 0, 0, 117.21471633]]) #21

In this case, np.linalg.det(M) returns -9504863423.43. I'm quite sure determinant is neither 0.0 nor -9504863423.43 because I calculated the same determinant with MATLAB and some online calculators as -3.8108e+019.

I also tried to perform LU decomposition and calculating the determinant with mpmath but these didn't work neither. How come these two cases may return different values? and any ideas how to calculate the determinant correctly?

I'm using Python(x,y) 2.7.6.1 on a 32-bit Win7 operating system.

ali_m
  • 71,714
  • 23
  • 223
  • 298
erjohn
  • 21
  • 4
  • 6
    Are you sure there is not a typo in the matrix you wrote? – Tom Cornebize Mar 26 '15 at 12:57
  • 1
    Sympy finds a determinant of -9504928268.73573 (I find the difference with numpy strange, but it is much closer to your result than -3.8108e19). – Tom Cornebize Mar 26 '15 at 12:58
  • 1
    Mathematica confirms the determinant to be -9.50486*10^9, and I also reproduced your numpy calculation. I think the full numpy matrix is fine. I'm not sure what's happening with the zero determinant, maybe some way you'rs setting the m's? – WakkaDojo Mar 26 '15 at 15:23
  • Also `scipy.sparse.linal.eigs` confirms a result similar to `-9504928268.73573` (k=20). But I'm not sure the if the algorithmic implementation differs. (py3.7-64bit, I only have this version) – Matteo Ragni Mar 16 '20 at 08:10

1 Answers1

0

The most precise and automated way will be to:

  1. Make sure/set matrix to use dtype='float32'
  2. Then do this:

    M =  # Your matrix
    (sign, logdet) = np.linalg.slogdet(M)
    determinant = np.exp(logdet)
    
Matteo Ragni
  • 2,837
  • 1
  • 20
  • 34
  • From a test it results that the determinant with type `np.float32` is `4.33165e+18`, while with `np.float64` is `-9504863423.42953`. Using `np.linalg.slogdet` returns the same results, loosing the sign. To recover the sign: `np.exp(logdet) * sign`. – Matteo Ragni Mar 16 '20 at 08:02