- 105 Views
- Uploaded on

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

Linear Systems – Gaussian Elimination

Linear Systems

- Solve Ax=b, where A is an nnmatrix andb is an n1 column vector
- We can also talk about non-square systems whereA is mn, b is m1, and x is n1
- 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

- If A is upper triangular, solve by back-substitution:

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.

Download Presentation

Connecting to Server..