4 ms·

I get that inverting a matrix can be computationally intensive, but what's the alternative? "What if you have to solve Ax = b for a lot of different b‘s? Surel

• askee 10y ago
  You would solve e.g. min_x* ( ||A x* - b|| (+\lamba regularization(x*) )
• jacobolus 10y ago
  For example, you factorize A into LU, so that you have LUx = b. Now you can let Ux = y, solve the system Ly = b for y, then solve Ux = y for x. Solving each triangular system is fast: O(n^2). If you need to do this with a different b, you can do it again, and you already have A factorized, which was the expensive part, i.e. O(n^3). (Depending on the type of problem, you might want to factorize the original matrix in a different way, e.g. QR or Cholesky or whatever.)
  • tclmeelmo 10y ago
    IIRC, the HP-15C used LU decomposition in its matrix math routines.
  • deleted 10y ago
    [deleted]
• brandmeyer 10y ago
  Sure, assuming that A is full-rank and small, then its perfectly reasonable to invert A with adjoint(A)/discriminant(A). But this is very expensive, and sensitive to roundoff for larger A. Lets say you are using QR factorization. Then R\b is easy to compute by back-substitution, because R is triangular. Q\b is also easy, because Q is orthogonal, and Q's inverse is its transpose, so you just calculate transpose(Q)*b. Same thing for LDL^T factorization, LU factorization, etc. Its easy because back-substitution with a triangular matrix is easy. In the case of LDL^T factorization, D is diagonal, and again D\b is easy (its a bunch of rank-1 division operations).
• neutronicus 10y ago
  In PDEs, we typically use something called a Krylov subspace iterative solver. These solvers construct a solution by adding up terms of the form c_n * A^n * y and truncating the series when the residuals are small enough. This is especially efficient in PDEs, where problem matrices are sparse so matrix-vector products are efficient (i.e. O(n)), and in practice only a few terms of the series are necessary.
• IngoBlechschmid 10y ago
  "How else would you solve Ax = b for x, if not by computing the inverse of A (assuming A is a 3x3 or 4x4)?" For instance, there a fast iterative solvers. They start with an initial guess for the solution (which can be an arbitrary vector) and successively generate more accurate approximations. These iterative solvers don't need access to the inverse of A. Furthermore, many of those solvers don't even need access to A: They only need access to a subroutine which, given a vector u, returns Au. Many matrices arising in applications (for instance one often reduces the problem of solving partial differential equations to the problem of solving a (very huge) linear system of equations) are sparse (have only few non-zero entries). For these it can make sense to write specialized code which specifically calculates Au. If you want to read more about this, then a good starting point is https://en.wikipedia.org/wiki/Kaczmarz_method https://en.wikipedia.org/wiki/Kaczmarz_method.
• donuteaters 10y ago
  You calculate the factorization (QR, LU, or similar), but you never calculate the inverse of the factorization: A x = b A => QR # QR decomposition Q R x = b R x = Q' b # Q is inverted by transposing x = backsubstitute(R, Q' b) Q was the transpose of it's inverse, and so you don't get any rounding error when you transpose it. R didn't need to be inverted because it was in a form where back substitution lets you apply the inverse of R to be without needing the inverse of R.
  • deleted 10y ago
    [deleted]
  • pramodliv1 10y ago
    MATLAB's linsolve uses your strategy! (http://in.mathworks.com/help/matlab/ref/linsolve.html http://in.mathworks.com/help/matlab/ref/linsolve.html)
    • ryanmonroe 10y ago
      Linear models in R are also solved by QR factorization.
• dfbrown 10y ago
  In LAPACK you would use dgetrs to solve for x directly with the LU factorization.
• JadeNB 10y ago
  I think that the idea is not that you don't take steps that, taken together, amount to inverting `A`, but rather that you don't just try to invert `A`. Thus, I think that the idea is to do something like "Compute `A = LU`, and then invert `L` and `U`"—which of course can be used to invert `A`, but which isn't exactly the same as "Compute `A^{-1}`" (using Cramer's rule, or some other absurdity). Incidentally, while solving `Ax = b` in general is definitely equivalent to inverting `A` (simply by taking `b` to be basis vectors, so that the solutions `x` are the columns of the inverse matrix in that basis), it is certainly possible to solve it in particular cases without inverting `A`; for example, the solution to `Ax = b`, where `b` is the first column of `A` in a particular basis, is the first vector `x` of that basis—no inversion required (nor even computation of `b`)!
• coherentpony 10y ago
  >Great. But don't you then use the factorization to determine the inverse of A (LAPACK example: call dgetrf to get the factorization, then call dgetri to compute the inverse using the factorization)? No. >How else would you solve Ax = b for x? Call dgetrs instead of dgetri. You use the factorisation directly (not the inverse) to get x from b. This is done by calling dgetrf to get the factorisation, then calling dgetrs to get the solution. You don't need to call dgetri. I'd do an example, but there are other replies to your comment that go through the mathematics.
• drauh 10y ago
  http://www.amazon.com/Applied-Numerical-Linear-Algebra-Demmel/dp/0898713897/ http://www.amazon.com/Applied-Numerical-Linear-Algebra-Demme... Numerical linear algebra has a long and vast history. You really do need specialized information that your undergrad linear algebra classes would not have imparted. For instance, the condition number is important: https://en.wikipedia.org/wiki/Condition_number https://en.wikipedia.org/wiki/Condition_number
• hchenji 10y ago
  Use a Pseudoinverse or a decomposition instead. A great book to read is Matrix Computations by Golub/van Loan.