Scientific Computing

1 / 36

# Scientific Computing - PowerPoint PPT Presentation

Scientific Computing. Linear Systems – Gaussian Elimination. Linear Systems. Linear Systems. Solve Ax=b, where A is an n  n matrix and b is an n  1 column vector We can also talk about non-square systems where A is m  n , b is m  1, and x is n  1

I am the owner, or an agent authorized to act on behalf of the owner, of the copyrighted work described.

## PowerPoint Slideshow about ' Scientific Computing' - alaura

Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author.While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server.

- - - - - - - - - - - - - - - - - - - - - - - - - - E N D - - - - - - - - - - - - - - - - - - - - - - - - - -
Presentation Transcript
Scientific Computing

Linear Systems – Gaussian Elimination

Linear Systems
• Solve Ax=b, where A is an nnmatrix andb is an n1 column vector
• We can also talk about non-square systems whereA is mn, b is m1, and x is n1
• Overdetermined if m>n:“more equations than unknowns”
• Underdetermined if n>m:“more unknowns than equations”
Singular Systems
• A is singular if some row islinear combination of other rows
• Singular systems can be underdetermined:or inconsistent:
Using Matrix Inverses
• To solve Ax = b it would be nice to use the inverse to A, that is, A-1
• However, it is usually not a good idea to compute x=A-1b
• Inefficient
• Prone to roundoff error
Gaussian Elimination
• Fundamental operations:
• Replace one equation with linear combinationof other equations
• Interchange two equations
• Multiply one equation by a scalar
• These are called elementary row operations.
• Do these operations again and again to reduce the system to a “trivial” system
Triangular Form
• Two special forms of matrices are especially nice for solving Ax=b:
• In both cases, successive substitution leads to a solution
Triangular Form
• A is lower triangular
Triangular Form
• Solve by forward substitution:
Triangular Form
• Solve by forward substitution:
Triangular Form
• Solve by forward substitution:

Etc

Triangular Form
• If A is upper triangular, solve by back-substitution:
Triangular Form
• Solve by back-substitution:

Etc

Gaussian Elimination Algorithm
• Do elementary row operations on the augmented system [A|b] to reduce the system to upper triangular form.
• Then, use back-substitution to find the answer.
Gaussian Elimination
• Example:
• Augmented Matrix form:
Gaussian Elimination
• Row Ops: What do we do to zero out first column under first pivot?
• Zero out below second pivot:
Gaussian Elimination
• Back-substitute
Matlab Implementation
• Task: Implement Gaussian Elimination (without pivoting) in a Matlab M-file.
• Notes
• Input = Coefficient matrix A, rhs b
• Output = solution vector x
Matlab Implementation
• Class Exercise: We will go through this code line by line, using the example in Pav section 3.3.2 to see how the code works.
Matlab Implementation
• Matlab:
• >> A=[2 1 1 3; 4 4 0 7; 6 5 4 17; 2 -1 0 7]
• A =
• 2 1 1 3
• 4 4 0 7
• 6 5 4 17
• 2 -1 0 7
• >> b = [7 11 31 15]\'
• b =
• 7
• 11
• 31
• 15

>> gauss_no_pivot(A,b)

ans =

1.5417

-1.4167

0.8333

1.5000

Matlab Implementation
• Class Exercise: How would we change the Matlab function gauss_no_pivot so we could see the result of each step of the row reduction?
Gaussian Elimination - Pivoting
• Consider this system:
• We immediately run into a problem:we cannot zero out below pivot, or back-substitute!
• More subtle problem:
Gaussian Elimination - Pivoting
• Conclusion: small diagonal elements are bad!
• Remedy: swap the row with the small diagonal element with a row below, this is called pivoting
Gaussian Elimination - Pivoting
• Our Example:
• Swap rows 1 and 2:
• Now continue:
Gaussian Elimination - Pivoting
• Two strategies for pivoting:
• Partial Pivoting
• Scaled Partial Pivoting
Matlab – Partial Pivoting
• Partial Pivoting:
• At step k, we are working on kth row, pivot = Akk . Search for largest Aik in kth column below (and including) Akk . Let p = index of row containing largest entry.
• If p ≠ k, swap rows p and k.
• Continue with Gaussian Elimination.
Matlab – Partial Pivoting
• Finding largest entry in a column:
• >> A
• A =
• 2 1 1 3
• 4 4 0 7
• 6 5 4 17
• 2 -1 0 7
• >> [r,m] = max(A(2:4,2))
• r =
• 5
• m =
• 2
• Why isn’t m = 3?
Matlab – Partial Pivoting
• Swapping rows m and k:
• BB =
• 2 1 1 3
• 4 4 0 7
• 6 5 4 17
• 2 -1 0 7
• >> BB([1 3],:) = BB([3 1], :)
• BB =
• 6 5 4 17
• 4 4 0 7
• 2 1 1 3
• 2 -1 0 7
Matlab – Partial Pivoting
• Code change?
• All that is needed is in main loop (going from rows k ->n-1) we add
• % Find max value (M) and index (m) of max entry below AAkk
• [M,m] = max(abs(AA(k:n,k)));
• m = m + (k-1); % row offset
• % swap rows, if needed
• if m ~= k, AA([k m],:) = AA([m k],:);
• end
Matlab – Partial Pivoting
• Class Exercise
• Review example in Moler Section 2.6
Scaled Partial Pivoting

Pav, Section 3.2

Practice
• Class Exercise: We will work through the example in Pav, Section 3.2.2
Practice
• Class Exercise: Do Pav Ex 3.7 for partial pivoting (“naïve Gaussian Elimination”).
• Do Pav Ex 3.7 for scaled partial pivoting.