Gaussian Elimination

• The procedure just described has a major problem.
• While it may be satisfactory for hand computations with small systems, it is inadequate for a large system.
• Observe that the transformed coefficients can become very large (getting bigger!) as we convert to a triangular system.
• The method that is called Gaussian elimination avoids this by subtracting $a_{i1}/a_{11}$ times the first equation from the $i^{th}$ equation to make the transformed numbers in the first column equal to zero.
• We do similarly for the rest of the columns.
• Observe that zeros may be created in the diagonal positions even if they are not present in the original matrix of coefficients.
• A useful strategy to avoid (if possible) such zero divisors is to rearrange the equations so as to put the coefficient of largest magnitude on the diagonal at each step.
• This is called pivoting.

• Repeat the example of the previous section,
  $$\left[ \begin{array}{rrrr} 4 & -2 & 1 &15 \\ -3 & -1 & 4 &8 \\ 1 & -1 & 3 &13 \\ \end{array} \right],$$
  $$\begin{array}{r} \\ R_2-(-3/4)R_1 \rightarrow \\ R_3-(1/4)R_1 \rightarrow \\ \end{array}$$
  $$\left[ \begin{array}{rrrr} 4 & -2 & 1 &15 \\ 0 & -2.5 & 4.75 &19.25 \\ 0 & -0.5 & 2.75 & 9.25 \\ \end{array} \right],$$
  $$\begin{array}{r} \\ \\ R_3-(-0.5/-2.5)R_2 \rightarrow \\ \end{array}$$
  $$\left[ \begin{array}{rrrr} 4 & -2 & 1 &15 \\ 0 & -2.5 & 4.75 & 19.25 \\ 0 & 0 & 1.8 & 5.40 \\ \end{array} \right]$$

• The method we have just illustrated is called Gaussian elimination.
• In this example, no pivoting was required to make the largest coefficients be on the diagonal.
• Back-substitution, gives us $x_3= 3, x_2=-2, x_1=2$

Example m-file: Show steps in Gauss elimination and back substitution. No pivoting. (http://siber.cankaya.edu.tr/ozdogan/NumericalComputations/mfiles/chapter2/GEshow.m GEshow.m)

• We shall obtain answers that are just close approximations to the exact answer because of round-off error.
• When there are many equations, the effects of round-off (the term is applied to the error due to chopping as well as when rounding is used) may cause large effects.
• In certain cases, the coefficients are such that the results are particularly sensitive to round off; such systems are called ill-conditioned.

• if we had stored the ratio of coefficients in place of zero (we show these in parentheses), our final form would have been
  $$\left[ \begin{array}{rrrr} 4 & -2 & 1 &15 \\ (-0.75) & -2... ...9.25 \\ (0.25) & (0.20) & 1.8 & 5.40 \\ \end{array} \right]$$

• The original matrix can be written as the product:

  
  

  $$A=$$
  $$\underbrace{\left[ \begin{array}{rrr} 1 & 0 & 0 \\ -0.75 & 1 & 0 \\ 0.25 & 0.20 & 1 \\ \end{array} \right]}_L*$$
  $$\underbrace{\left[ \begin{array}{rrr} 4 & -2 & 1 \\ 0 & -2.5 & 4.75 \\ 0 & 0 & 1.8 \\ \end{array} \right]}_U$$

• This procedure is called a $LU$-decomposition of $A$.
• In this case,
  $$A =L *U$$
  where $L$ is lower-triangular and $U$ is upper-triangular.

• The determinant of two matrices, $B \times C$, is the product of each of the determinants, for this example we have
  $$det(L*U) = det(L) * det(U)= det(U)$$
  Because $L$ is triangular and has only ones on its diagonal so that $det(L)=1$.
• Thus, for our example, we have
  $$det(A) = det(U)=(4)*(-2.5)*(1.8)=-18$$
  since $U$ is upper-triangular and its determinant is just the product of the diagonal elements.
• When there are row interchanges
  $$det(A) = (-1)^m*u_{11}*\ldots*u_{nn}$$
  where the exponent m represents the number of row interchanges.

• Example. Solve the following system of equations using Gaussian elimination.
  $$\begin{array}{rrrr} &2x_2& &+x_4=0\\ 2x_1 &+2x_2 &+3x_3 &... ..._2 & & x_4 =-7 \\ 6x_1 &+x_2 &-6x_3 &-5x_4 =6 \\ \end{array}$$

• In addition, obtain the determinant of the coefficient matrix and the $LU$ decomposition of this matrix.

1. The augmented coefficient matrix is
   $$\left[ \begin{array}{rrrrr} 0 & 2 & 0 & 1 &0 \\ 2 & 2 & 3... ...& -3 & 0 & 1 &-7 \\ 6 & 1 &-6 &-5 &6 \\ \end{array} \right]$$

2. We cannot permit a zero in the $a_{11}$ position because that element is the pivot in reducing the first column.
3. We could interchange the first row with any of the other rows to avoid a zero divisor, but interchanging the first and fourth rows is our best choice. This gives
   $$\left[ \begin{array}{rrrrr} 6 & 1 &-6 &-5 &6 \\ 2 & 2 & 3... ...& -3 & 0 & 1 &-7 \\ 0 & 2 & 0 & 1 &0 \\ \end{array} \right]$$
   $$\left[ \begin{array}{rrrrr} 6 & 1 &-6 &-5 &6 \\ 0 &1.6667... ... & 4 &4.3333&-11 \\ 0 & 2 & 0 & 1 &0 \\ \end{array} \right]$$

4

We again interchange before reducing the second column, not because we have a zero divisor, but because we want to preserve accuracy. Interchanging the second and third rows puts the element of largest magnitude on the diagonal.

$$\left[ \begin{array}{rrrrr} 6 & 1 &-6 &-5 &6 \\ 0 &-3.666... ...67& 5 &3.6667&-4 \\ 0 & 2 & 0 & 1 &0 \\ \end{array} \right]$$

5

Now we reduce in the second column

$$\left[ \begin{array}{rrrrr} 6 & 1 &-6 &-5 &6 \\ 0 &-3.666... ...1\\ 0 & 0 & 2.1818 & 3.3636 &-5.9999 \\ \end{array} \right]$$

6

No interchange is indicated in the third column. Reducing, we get

$$\left[ \begin{array}{rrrrr} 6 & 1 &-6 &-5 &6 \\ 0 &-3.666... ...9.0001\\ 0 & 0 & 0 & 1.5600 &-3.1199 \\ \end{array} \right]$$

7

Back-substitution gives

$$x_1=-0.50000, x_2=1.0000, x_3=0.33325, x_4=-1.9999.$$

• The correct (exact) answers are $x_1=-1/2,x_2=1,x_3=1/3,x_4=-2$.
• In this calculation we have carried five significant figures and rounded each calculation.
• Even so, we do not have five-digit accuracy in the answers.
• The discrepancy is due to round off.
• Example m-files: (http://siber.cankaya.edu.tr/ozdogan/NumericalComputations/mfiles/chapter2/GEshow.m GEshow.m, http://siber.cankaya.edu.tr/ozdogan/NumericalComputations/mfiles/chapter2/GEPivShow.m GEPivShow.m)

  

• In this example, if we had replaced the zeros below the main diagonal with the ratio of coefficients at each step, the resulting augmented matrix would be
  $$\left[ \begin{array}{rrrrr} 6 & 1 &-6 &-5 &6 \\ (0.66667)... ... (-0.54545) & (0.32) & 1.5600 &-3.1199 \\ \end{array} \right]$$

• This gives a LU decomposition as
  $$\left[ \begin{array}{rrrrr} 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]$$

• It should be noted that the product of these matrices produces a permutation of the original matrix, call it $A'$, where
  $$A'=\left[ \begin{array}{rrrr} 6 & 1 &-6 &-5 \\ 4 & -3 & 0... ... \\ 2 & 2 & 3 & 2 \\ 0 & 2 & 0 & 1 \\ \end{array} \right]$$

• The determinant of the original matrix of coefficients can be easily computed according to the formula
  $$det(A)=(-1)^2*(6)*(-3.6667)*(6.8182)*(1.5600)=-234.0028$$
  which is close to the exact solution: -234.
• The exponent 2 is required, because there were two row interchanges in solving this system.
• To summarize
  1. The solution to the four equations
  2. The determinant of the coefficient matrix
  3. A $LU$ decomposition of the matrix, $A'$, which is just the original matrix, $A$, after we have interchanged its rows.
• ``These'' are readily obtained after solving the system by Gaussian elimination method.

Example m-file: LU factorization without pivoting. (http://siber.cankaya.edu.tr/ozdogan/NumericalComputations/mfiles/chapter2/luNopiv.m luNopiv.m)
Example m-file: LU factorization with partial pivoting. (http://siber.cankaya.edu.tr/ozdogan/NumericalComputations/mfiles/chapter2/luPiv.m luPiv.m)