Using the LU Matrix for Multiple Right-Hand Sides

• Many physical situations are modelled with a large set of linear equations.
• The equations will depend on the geometry and certain external factors that will determine the right-hand sides.
• If we want the solution for many different values of these right-hand sides,
• it is inefficient to solve the system from the start with each one of the right-hand-side values.
• Using the $LU$ equivalent of the coefficient matrix is preferred.

• Suppose we have solved the system $Ax=b$ by Gaussian elimination.
• We now know the $LU$ equivalent of $A$:
  $A=L*U$
• Consider now that we want to solve $Ax=b$ with some new $b$-vector.
• We can write
  
  $$Ax=b$$
  
  
  
  
  $$LUx=b$$
  
  
  
  
  $$Ly=b$$
  
  

• The product of $U$ and $x$ is a vector, call it $y$.
• Now, we can solve for $y$ from $Ly = b$.
• This is readily done because $L$ is lower-triangular and we get $y$ by forward-substitution.
• Call the solution $y = b'$.
• Going back to the original $LUx = b$, we see that, from $Ux=y=b'$, we can get $x$ from $Ux = b'$.
• Which is again readily done by back-substitution ($U$ is upper-triangular).

• i.e., Solve $Ax=b$, where we already have its $L$ and $U$ matrices:

  
  

  $$\left[ \begin{array}{rrrr} 1&0&0&0\\ 0.66667 & 1 & 0 & 0 ... ... & 1 &0 \\ 0.0 & -0.54545 & 0.32 & 1 \\ \end{array}\right]*$$
  
  
  
  
  $$\left[ \begin{array}{rrrr} 6 & 1 &-6 &-5 \\ 0&-3.6667 & ... ... & 6.8182 &5.6364\\ 0& 0 &0 & 1.5600 \\ \end{array} \right]$$
  
  

• Suppose that the $b$-vector is $[6~ -7~ -2~ 0]^T$.
• We first get $y=Ux$ from $Ly = b$ by forward substitution:
  
  $$y = [6~ -11~-9~ -3.12]^T$$
  
  

• and use it to compute $x$ from $Ux = y$:
  
  $$x= [ -0.5~ 1~ 0.3333~ -2]^T.$$
  
  

• Now, if we want the solution with a different $b$-vector;
  
  $$bb = [1~ 4~ -3~ 1]^T$$
  
  

• we just do $Ly = bb$ to get
  
  $$y = [1~ 3.3333~-1.8182~3.4]^T$$
  
  

• and then use this $y$ in $Ux = y$ to find the new $x$:
  
  $$x=[0.0128~ -0.5897~ -2.0684~ 2.1795]^T$$