scientific computing
Download
Skip this Video
Download Presentation
Scientific Computing

Loading in 2 Seconds...

play fullscreen
1 / 36

Scientific Computing - PowerPoint PPT Presentation


  • 105 Views
  • Uploaded on

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

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

PowerPoint Slideshow about ' Scientific Computing' - alaura


An Image/Link below is provided (as is) to download presentation

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
Scientific Computing

Linear Systems – Gaussian Elimination

linear systems1
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
Singular Systems
  • A is singular if some row islinear combination of other rows
  • Singular systems can be underdetermined:or inconsistent:
using matrix inverses
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
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
Triangular Form
  • Two special forms of matrices are especially nice for solving Ax=b:
  • In both cases, successive substitution leads to a solution
triangular form1
Triangular Form
  • A is lower triangular
triangular form2
Triangular Form
  • Solve by forward substitution:
triangular form3
Triangular Form
  • Solve by forward substitution:
triangular form4
Triangular Form
  • Solve by forward substitution:

Etc

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

Etc

gaussian elimination algorithm
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 elimination1
Gaussian Elimination
  • Example:
  • Augmented Matrix form:
gaussian elimination2
Gaussian Elimination
  • Row Ops: What do we do to zero out first column under first pivot?
  • Zero out below second pivot:
gaussian elimination3
Gaussian Elimination
  • Back-substitute
matlab implementation
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 implementation1
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 implementation2
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 implementation3
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
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 pivoting1
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 pivoting2
Gaussian Elimination - Pivoting
  • Our Example:
  • Swap rows 1 and 2:
  • Now continue:
gaussian elimination pivoting3
Gaussian Elimination - Pivoting
  • Two strategies for pivoting:
    • Partial Pivoting
    • Scaled Partial Pivoting
matlab 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 pivoting1
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 pivoting2
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 pivoting3
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 pivoting4
Matlab – Partial Pivoting
  • Class Exercise
  • Review example in Moler Section 2.6
scaled partial pivoting
Scaled Partial Pivoting

Pav, Section 3.2

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