1 / 31

Engineering Computation

Engineering Computation. Lecture 7. Errors in Solutions to Systems of Linear Equations. System of Equations Errors in Solutions to Systems of Linear Equations Objective: Solve [A]{x} = {b} Problem:

Download Presentation

Engineering Computation

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. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript


  1. EngineeringComputation Lecture 7 E. T. S. I. Caminos, Canales y Puertos

  2. Errors in Solutions to Systems of Linear Equations • System of Equations • Errors in Solutions to Systems of Linear Equations • Objective:Solve [A]{x} = {b} • Problem: • Round-off errors may accumulate and even be exaggerated by the solution procedure. Errors are often exaggerated if the system is ill-conditioned • Possible remedies to minimize this effect: • 1. Partial or complete pivoting • 2. Work in double precision • 3. Transform the problem into an equivalent system of linear equations by scaling or equilibrating E. T. S. I. Caminos, Canales y Puertos

  3. Errors in Solutions to Systems of Linear Equations • Ill-conditioning • A system of equations is singular if det|A| = 0 • If a system of equations is nearly singular it is ill-conditioned. • Systems which are ill-conditioned are extremely sensitive to small changes in coefficients of [A] and {b}. These systems are inherently sensitive to round-off errors. • Question: • Can we develop a means for detecting these situations? E. T. S. I. Caminos, Canales y Puertos

  4. a11x1+ a12x2 = b1 x2 b2/a22 b2/a21 b1/a12 x1 b1/a11 a21x1+ a22x2 = b2 Errors in Solutions to Systems of Linear Equations Ill-conditioning of [A]{x} = {b}: Consider the graphical interpretation for a 2-equation system: We can plot the two linear equations on a graph of x1 vs. x2. E. T. S. I. Caminos, Canales y Puertos

  5. Errors in Solutions to Systems of Linear Equations Ill-conditioning of [A]{x} = {b}: Consider the graphical interpretation for a 2-equation system: We can plot the two linear equations on a graph of x1 vs. x2. x1 x1 x2 x2 Uncertainty in x2 Uncertainty in x2 Well-conditioned Ill-conditioned E. T. S. I. Caminos, Canales y Puertos

  6. Errors in Solutions to Systems of Linear Equations Ways to detect ill-conditioning: 1. Calculate {x}, make small change in [A] or {b} and determine the change in the solution {x}. 2. After forward elimination, examine diagonal of upper triangular matrix. If aii << ajj, i.e. there is a relatively small value on diagonal, then this may indicate ill-conditioning. 3. Compare {x}single with {x}double 4. Estimate the "condition number" for A. Substituting the calculated {x} into [A]{x} and checking this against {b} will not always work!!! E. T. S. I. Caminos, Canales y Puertos

  7. Errors in Solutions to Systems of Linear Equations • Ways to detect ill-conditioning: • If det|A| = 0 the matrix is singular • ==> the determinant may be an indicator of conditioning • If det|A| is near zero is the matrix ill-conditioned? • Consider: • After scaling: • ==> det|A| will provide an estimate of conditioning if it is normalized by the "magnitude" of the matrix. E. T. S. I. Caminos, Canales y Puertos

  8. Norms Norms and the Condition Number We need a quantitative measure of ill-conditioning. This measure will then directly reflect the possible magnitude of round off effects. To do this we need to understand norms: Norm:Scalar measure of the magnitude of a matrix or vector ("how big" a vector is). Not to be confused with the dimension of a matrix. E. T. S. I. Caminos, Canales y Puertos

  9. Vector Norms Vector Norms: Scalar measure of the magnitude of a vector Here are some vector norms for n x 1 vectors {x} with typical elements xi. Each is in the general form of a p norm defined by the general relationship: 1. Sum of the magnitudes: 2. Magnitude of largest element: (infinity norm) 3. Length or Euclidean norm: E. T. S. I. Caminos, Canales y Puertos

  10. Norms • Vector Norms • Required Properties of vector norm: • 1. ||x||  0 and ||x|| = 0 if and only if [x]=0 • 2 ||kx|| = k ||x|| where k is any positive scalar • 3. ||x+y||  ||x|| + ||y|| Triangle Inequality • For the Euclidean vector norm we also have • 4. ||x•y||  ||x|| ||y|| • because the dot product or inner product property satisfies: • ||xy|| = ||x||•||y|| |cos()|  ||x|| • ||y||. E. T. S. I. Caminos, Canales y Puertos

  11. Matrix Norms Matrix Norms: Scalar measure of the magnitude of a matrix. Matrix norms corresponding to vector norms above are defined by the general relationship: 1. Largest column sum: (column sum norm) 2. Largest row sum: (row sum norm) (infinity norm) E. T. S. I. Caminos, Canales y Puertos

  12. Matrix norms • 3. Spectral norm: ||A|| 2 = (µmax)1/2 • where µmax is the largest eigenvalue of [A]T[A] • If [A] is symmetric, (µmax)1/2 = max , is the largest eigenvalue of [A]. • (Note: this is not the same as the Euclidean or Frobenius norm, seldom used: E. T. S. I. Caminos, Canales y Puertos

  13. Matrix norms • MatrixNorms • For matrix norms to be useful we require that • 0. || Ax ||  || A || ||x || • General properties of any matrix norm: • 1. || A ||  0 and || A || = 0 iff [A] = 0 • 2. || k A || = k || A || where k is any positive scalar • 3. || A + B ||  || A || + || B || "Triangle Inequality" • 4. || A B ||  || A || || B || • Why are norms important? • Norms permit us to express the accuracy of the solution {x} in terms of ||x|| • Norms allow us to bound the magnitude of the product [ A ] {x} and the associated errors. E. T. S. I. Caminos, Canales y Puertos

  14. Error Analysis • Forward and backward error analysis can estimate the effect of truncation and roundoff errors on the precision of a result. The two approaches are alternative views: • Forward (a priori) error analysis tries to trace the accumulation of error through each process of the algorithm, comparing the calculated and exact values at every stage. • Backward (a posteriori) error analysis views the final solution as the exact solution to a perturbed problem. One can consider how different the perturbed problem is from the original problem. • Here we use the condition number of a matrix [A] to specify the amount by which relative errors in [A] and/or {b} due to input, truncation, and rounding can be amplified by the linear system in the computation of {x}. E. T. S. I. Caminos, Canales y Puertos

  15. Error Analysis • Backward Error Analysis of [A]{x} = {b} for errors in {b} • Suppose the coefficients {b} are not precisely represented. What might be the effect on the calculated value for {x + dx}? • Lemma: [A]{x} = {b} yields ||A|| ||x||  ||b|| or • Now an error in {b} yields a corresponding error in {x}: • [A ]{x + dx} = {b + db} • [A]{x} + [A]{ dx} = {b} + {db} • Subtracting [A]{x} = {b} yields: • [A]{dx} = {db} ––> {dx} = [A]-1{db} E. T. S. I. Caminos, Canales y Puertos

  16. Error Analysis Backward Error Analysis of [A]{x} = {b} for errors in {b} Taking norms we have: And using the lemma: we then have : Define the condition number as k = cond [A]  ||A-1|| ||A||  1 If k 1 or k is small, the system is well-conditioned If k >> 1, system is ill conditioned. 1 = || I || = || A-1A ||  || A-1 || || A || = k = Cond(A) E. T. S. I. Caminos, Canales y Puertos

  17. Error Analysis Backward Error Analysis of [A]{x} = {b} for errors in [A] If the coefficients in [A] are not precisely represented, what might be effect on the calculated value of {x+ dx}? [A + dA ]{x + dx} = {b} [A]{x} + [A]{ dx} + [dA]{x+dx} = {b} Subtracting [A]{x} = {b} yields: [A]{ dx} = – [dA]{x+dx} or {dx} = – [A]-1 [dA] {x+dx} Taking norms and multiplying by || A || / || A || yields : E. T. S. I. Caminos, Canales y Puertos

  18. Error Analysis Estimate of Loss of Significance: Consider the possible impact of errors [dA] on the precision of {x}. • implies that if • Or, taking log of both sides: s > p - log10() • log10() is the loss in decimal precision; i.e., we start with p decimal figures and end-up with s decimal figures. • It is not always necessary to find [A]-1 to estimate k = cond[A]. Instead, use an estimate based upon iteration of inverse matrix using LU decomposition. E. T. S. I. Caminos, Canales y Puertos

  19. Iterative Solution Methods • Impetus for Iterative Schemes: • 1. May be more rapid if coefficient matrix is "sparse" • 2. May be more economical with respect to memory • 3. May also be applied to solve nonlinear systems • Disadvantages: • 1. May not converge or may converge slowly • 2. Not appropriate for all systems Error bounds apply to solutions obtained by direct and iterative methods because they address the specification of [dA] and {db}. E. T. S. I. Caminos, Canales y Puertos

  20. Iterative Solution Methods Basic Mechanics: Starting with: a11x1 + a12x2 + a13x3 + ... + a1nxn = b1 a21x1 + a22x2 + a23x3 + ... + a2nxn = b2 a31x1 + a32x2 + a33x3 + ... + a3nxn = b3 : : an1x1 + an2x2 + an3x3 + ... + annxn = bn Solve each equation for one variable: x1 = [b1 – (a12x2 + a13x3 + ... + a1nxn )} / a11 x2 = [b2 – (a21x1 + a23x3 + ... + a2nxn )} / a22 x3 = [b3 – (a31x1 + a32x2 + ... + a3nxn )} / a33 : xn = [bn – (an1x2 + an2x3 + ... + an,n-1xn-1 )} / ann E. T. S. I. Caminos, Canales y Puertos

  21. Iterative Solution Methods • Start with initial estimate of {x}0. • Substitute into the right-hand side of all the equations. • Generate new approximation {x}1. • This is a multivariate one-point iteration: {x}j+1 = {g({x}j)} • Repeat process until the maximum number of iterations is reached or until: || xj+1 – xj ||  d + e || xj+1 || E. T. S. I. Caminos, Canales y Puertos

  22. Convergence • To solve [A]{x} = {b} • Separate [A] into: [A] = [Lo] + [D] + [Uo] • [D] = diagonal (aii) • [Lo] = lower triangular with 0's on diagonal • [Uo] = upper triangular with 0's on diagonal • Rewrite system: • [A]{x} = ( [Lo] + [D] + [Uo] ){x} = {b} • [D]{x} + ( [Lo] + [Uo] ){x} = {b} • Iterate: • [D]{x}j+1 = {b} – ( [Lo]+[Uo] ) {x}j • {x}j+1 = [D]-1{b} – [D]-1 ( [Lo]+[Uo] ) {x}j • Iterations converge if: • || [D]-1 ( [Lo] + [Uo] ) || < 1 • (sufficient if equations are diagonally dominant) E. T. S. I. Caminos, Canales y Puertos

  23. Iterative Solution Methods – the Jacobi Method E. T. S. I. Caminos, Canales y Puertos

  24. Iterative Solution Methods -- Gauss-Seidel In most cases using the newest values within the right-hand side equations will provide better estimates of the next value. If this is done, then we are using the Gauss-Seidel Method: ( [Lo]+[D] ){x}j+1 = {b} – [Uo] {x}j or explicitly: If this is done, then we are using the Gauss-Seidel Method E. T. S. I. Caminos, Canales y Puertos

  25. Iterative Solution Methods -- Gauss-Seidel If either method is going to converge, Gauss-Seidel will converge faster than Jacobi. Why use Jacobi at all? Because you can separate the n-equations into n independent tasks, it is very well suited computers with parallel processors. E. T. S. I. Caminos, Canales y Puertos

  26. Convergence of Iterative Solution Methods • Rewrite given system: [A]{x} = { [B] + [E] } {x} = {b} • where [B] is diagonal, or triangular so we can solve [B]{y} = {g} quickly. Thus, • [B] {x}j+1= {b}– [E] {x}j • which is effectively: {x}j+1 = [B]-1 ({b} – [E] {x}j ) • True solution {x}c satisfies: {x}c = [B]-1 ({b} – [E] {x}c) • Subtracting yields: {x}c – {x}j+1= – [B]-1 [E] [{x}c – {x}j] • So ||{x}c – {x}j+1 ||  ||[B]-1 [E]|| ||{x}c – {x}j || • Iterations converge linearly if || [B]-1 [E] || < 1=> || ([D] + [Lo])-1 [Uo] || < 1 For Gauss-Seidel • => || [D] -1 ([Lo] + [Uo]) || < 1 For Jacobi E. T. S. I. Caminos, Canales y Puertos

  27. Convergence of Iterative Solution Methods • Iterative methods will not converge for all systems of equations, nor for all possible rearrangements. • If the system is diagonally dominant, • i.e., | aii | > | aij | where i  j then with all < 1.0, i.e., small slopes. E. T. S. I. Caminos, Canales y Puertos

  28. Convergence of Iterative Solution Methods A sufficient condition for convergence exists: • Notes: • 1. If the above does not hold, still may converge. • 2. This looks similar to infinity norm of [A] E. T. S. I. Caminos, Canales y Puertos

  29. Improving Rate of Convergence of G-S Iteration • Relaxation Schemes: • where 0.0 < l 2.0 • (Usually the value of l is close to 1) • Underrelaxation ( 0.0 < l< 1.0 ) • More weight is placed on the previous value. Often used to: • - make non-convergent system convergent or • - to expedite convergence by damping out oscillations. • Overrelaxation ( 1.0 < l 2.0 ) • More weight is placed on the new value. Assumes that the new value is heading in right direction, and hence pushes new value close to true solution. • The choice of l is highly problem-dependent and is empirical, so relaxation is usually only used for often repeated calculations of a particular class. E. T. S. I. Caminos, Canales y Puertos

  30. Why Iterative Solutions? • We often need to solve [A]{x} = {b} where n = 1000's • • Description of a building or airframe, • • Finite-Difference approximations to PDE's. Most of A's elements will be zero; a finite-difference approximation to Laplace's equation will have five aij0 in each row of A. • Direct method (Gaussian elimination) • • Requires n3/3 flops (say n = 5000; n3/3 = 4 x 1010 flops) • • Fills in many of n2-5n zero elements of A • Iterative methods (Jacobi or Gauss-Seidel) • • Never store [A] (say n = 5000; [A] would need 4n2 = 100 Mb) • • Only need to compute [A-B] {x}; and to solve [B]{xt+1} = {b} E. T. S. I. Caminos, Canales y Puertos

  31. Why Iterative Solutions? • • Effort: • Suppose [B] is diagonal, • solving [B] {v} = {b} n flops • Computing [A-B] x 4n flops • For m iterations 5mn flops • For n = m = 5000, 5mn = 1.25x108 • At worst O(n2). E. T. S. I. Caminos, Canales y Puertos

More Related