1. PF Contest - Win "Conquering the Physics GRE" book! Click Here to Enter
    Dismiss Notice
Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Matrix Inversion Problem

  1. Jun 11, 2012 #1
    1. The problem statement, all variables and given/known data

    Hi, I'm trying to calculate the inverse of a really big matrix (4096x4096) using matlab. The inversion process takes too long on my pc, so i want to find a faster way.
    The matrix has a strange property that i think could be useful to solve the problem:
    splitting it in 64 blocks (8 rows and 8 columns of 256x256 submatrices), each submatrice is diagonal. Inverting a diagonal matrix is really simple, so i am asking you if this property can be used to calculate the inversion in a faster way.

    2. Relevant equations

    3. The attempt at a solution

    i found something about block diagonal matrix inversion, but my matrix is a bit different from a block diagonal matrix, because i have no zero submatrices... each block is diagonal...
    i'm trying to use the matrix inversion in block form but i feel that there must be a smarter solution...
  2. jcsd
  3. Jun 11, 2012 #2

    Ray Vickson

    User Avatar
    Science Advisor
    Homework Helper

    You are correct: the special structure allow you to reduce the computations from many billions of multiplications and additions, to a few hundred thousand. It gets complicated (but straightforward in principle), so I will just illustrate it for two cases: (I) the 2n x 2n case (4 nxn diagonal blocks); and (II) the 3n x 3n case (9 nxn diagonal blocks).
    (I) [tex]M = \pmatrix{A&B\\D&E}, \text{ where } A = \text{diag}(a_1, a_2,\ldots,a_n), \text{ etc.}[/tex] We want to solve the equations [itex] MZ = R,[/itex] where
    [tex]R = \pmatrix{u\\v}, u, \, v = \text{given n-dimensional vectors}[/tex] and
    [tex] Z = \pmatrix{x\\y}, \; x,\, y = \text{unknown n-dimensional vectors.}[/tex]
    These become
    [tex] Ax + By = u \longrightarrow x = A^{-1}u - A^{-1}B y\\
    Dx+Ey = v \longrightarrow D(A^{-1}u - A^{-1}B y) + Ey = v,\\
    \text{so } E'y = v', \text{where } E' = E - D A^{-1}B, \; v' = v - D A^{-1} u.
    So, we can get y, then can get x from a previous formula.

    Note: this is simplified by using the following simple facts: (1) The inverse of a diagonal matrix is diagonal (if it exists), that is: [tex] (\text{diag}(a_1,\ldots,a_n))^{-1} = \text{diag}(1/a_1, \ldots, 1/a_n).[/tex] (2) The product of two diagonal matrices is diagonal; that is:
    [tex] \text{diag}(a_1,\ldots,a_n) \cdot \text{diag}(b_1,\ldots,b_n)
    = \text{diag}(a_1 b_1, \ldots, a_n b_n).[/tex].
    (3) The product of a scalar and a diagonal matrix is diagonal; that is
    [tex] c\: \text{diag}(a_1,\ldots,a_n) = \text{diag}(ca_1,\ldots, ca_n),[/tex] for any constant c.
    (4) The sum of diagonal matrices is diagonal; that is
    [tex] \text{diag}(a_1,\ldots,a_n) + \text{diag}(b_1,\ldots,n_n)
    = \text{diag}(a_1 + b_1, \ldots, a_n + b_n). [/tex]

    Therefore, the matrix E' above is diagonal, with (i,i)-element equal to
    [tex] E'(i,i) \equiv e'_i = e_i - d_i b_i/a_i, [/tex] so no actual matrix multiplications are needed. Further, solving E'y = v' for y and Ax = u-By for x is easy.

    (II) Solve
    [tex] \pmatrix{A&B&C\\D&E&F\\G&H&K} \pmatrix{x\\y\\z} = \pmatrix{u\\v\\w},[/tex]
    where each of the listed sub-matrices are nxn diagonal.
    We can write the first two equations as
    [tex] Ax + By = U, Dx + Ey = V, [/tex] where [itex] U = u-Cz,\; V = v - Fz,[/itex] so from above we have
    [tex] E' y = V' = V - DA^{-1}U = v - DA^{-1}u + (DA^{-1}C -F)z,[/tex]
    which we can write as [itex] E' y = v' + F' z.[/itex]
    We can substitute this into the equation for x to get [itex] Ax = u' + C'z, [/itex] where we get u' and C' from messy but straightforward algebra. So now we have x and y solved in terms of z, and can substitute these expressions into the third equation
    [tex] Gx + Hy + Kz = w,[/tex]
    to get an equation of the form [itex] K' z = w'.[/itex] All the matrices E', F', C', K', etc., are diagonal, so all the computations are very easy.

    You can extend the method to 64 blocks in the same way; it just needs lots of labels to keep track of the different matrices throughout. All of the matrices encountered are diagonal, so it is enough to store their diagonals, and no acutal matrix multiplications or inversions are needed.

    Last edited: Jun 11, 2012
  4. Jun 11, 2012 #3
    Very nice post, Ray!!
  5. Jun 12, 2012 #4
    thank you very much for your help!
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Similar Threads - Matrix Inversion Problem Date
Finding matrix and the inverse is given May 17, 2017
Find inverse matrix using determinants and adjoints Jun 21, 2016
Inverse of Matrix Mar 4, 2016
Inverse matrix word problem, matrix arithmetic Sep 28, 2014
Matrix inverse problem Feb 27, 2006