# Necessity/Advantage of LU Decomposition over Gaussian Elimination

I am reading the book "Introduction to Linear Algebra" by Gilbert Strang and couldn't help wondering the advantages of LU decomposition over Gaussian Elimination!

For a system of linear equations in the form Ax = b, one of the methods to solve the unknowns is Gaussian Elimination, where you form a upper triangular matrix U by forward elimination and then figure out the unknowns by backward substitution. This serves the purpose of solving a system of linear equations. What was necessity for LU Decomposition, i.e. after finding U by forward elimination, why do we go about finding L (the lower triangular matrix) when you already had U and could have done a backward substitution?

In many engineering applications, when you solve $Ax = b$, the matrix $A \in \mathbb{R}^{N \times N}$ remains unchanged, while the right hand side vector $b$ keeps changing.

A typical example is when you are solving a partial differential equation for different forcing functions. For these different forcing functions, the meshing is usually kept the same. The matrix $A$ only depends on the mesh parameters and hence remains unchanged for the different forcing functions. However, the right hand side vector $b$ changes for each of the forcing function.

Another example is when you are solving a time dependent problem, where the unknowns evolve with time. In this case again, if the time stepping is constant across different time instants, then again the matrix $A$ remains unchanged and the only the right hand side vector $b$ changes at each time step.

The key idea behind solving using the $LU$ factorization (for that matter any factorization) is to decouple the factorization phase (usually computationally expensive) from the actual solving phase. The factorization phase only needs the matrix $A$, while the actual solving phase makes use of the factored form of $A$ and the right hand side to solve the linear system. Hence, once we have the factorization, we can make use of the factored form of $A$, to solve for different right hand sides at a relatively moderate computational cost.

The cost of factorizing the matrix $A$ into $LU$ is $\mathcal{O}(N^3)$. Once you have this factorization, the cost of solving i.e. the cost of solving $LUx = b$ is just $\mathcal{O}(N^2)$, since the cost of solving a triangular system scales as $\mathcal{O}(N^2)$.

(Note that to solve $LUx = b$, you first solve $Ly = b$ and then $Ux = y$. Solving $Ly = b$ and $Ux=y$ costs $\mathcal{O}(N^2).$)

Hence, if you have '$r$' right hand side vectors $\{b_1,b_2,\ldots,b_r\}$, once you have the $LU$ factorization of the matrix $A$, the total cost to solve $$Ax_1 = b_1, Ax_2 = b_2 , \ldots, Ax_r = b_r$$ scales as $\mathcal{O}(N^3 + rN^2)$.

On the other hand, if you do Gauss elimination separately for each right hand side vector $b_j$, then the total cost scales as $\mathcal{O}(rN^3)$, since each Gauss elimination independently costs $\mathcal{O}(N^3)$.

However, typically when people say Gauss elimination, they usually refer to $LU$ decomposition and not to the method of solving each right hand side completely independently.

• Thanks Marvis! That made it clear! :) – Chethan Ravindranath Dec 28 '12 at 13:27

We have the set of $n$ systems $A[x_1, x_2, ..., x_n] = [b_1, b_2, ..., b_n]$. Why use LU-decomposition?

1) Why not just solve this set of $n$ systems simultaniusly using Gaussian Elimination?
2)Peforming LU-decomposition why do we stop at the step $L^{-1}A = U$ (from which we get $A = LU$)? Why not go further and get $A^{-1}A = I$ (suppose $A$ is invertible)? And then easily get the solutions $[x_1, x_2, ..., x_n] = A^{-1}[b_1, b_2, ..., b_n]$.

• Getting $A^{-1}$ is a quite bad idea, especially if $A$ is sparse (hint: $A^{-1}$ is not). In addition, a simple $solve()$ is faster than the two steps of first getting $A^{-1}$ by the use of $solve()$, and then an additional multiplication. – arc_lupus Apr 8 '17 at 14:17

The computational cost of solving $$x$$ for $$Ax = b$$ via Gaussian elimination or $$LU-$$decomposition is the same. Solving $$Ax = b$$ via Gaussian elimination with partial pivoting on the augmented matrix $$[A | b]$$ transforms $$[A|b]$$ into an upper-triangular system $$[U | b^{'}]$$. Thereafter, backward substitution is performed to determine $$x$$. Furthermore, in the $$LU-$$decomposition matrix, $$L$$ is nothing but "multipliers" obtained during the Gauss elimination process. So there is no additional cost of computing $$L$$.

But consider a scenario when the right-hand side $$Ax = b$$ keeps on changing and $$A$$ is fixed. That means you deal with the system of equations $$Ax = b_1, Ax=b_2, ...Ax=b_k$$. In this case, $$LU-$$decomposition works more efficiently than Gauss elimination as you do not need to solve all augmented matrices. You just have to decompose $$A$$ into $$LU$$ (or Gaussian elimination with cost ~ $$n^3$$) once.