Skip to content
Textbook Guides

Solution Manual for Fundamentals of Electromagnetics with Engineering Applications, 1st Edition

Solution Manual for Fundamentals of Electromagnetics with Engineering Applications, 1st Edition helps you navigate your textbook with ease, offering answers to every question.

Preview 16 of 56 pages

See the first 16 pages exactly as they appear in the document.

  1. Page 1
    Read page 1 as text

    solutions MA nu A l F o R n umerical Algorithms Justin s olomon by 9781482251944_SM_Cover.indd 1 28/09/15 12:38 pm

  2. Page 2
  3. Page 3
    Read page 3 as text

     1 1.1. We compute the gradients element-by-element: ∇ f = ( ∂ ∂x ( x 2 + y 2 ) , ∂ ∂y ( x 2 + y 2 ) ) = (2 x, 2 y ) ∇ g = ( ∂ ∂x √ x 2 + y 2 , ∂ ∂y √ x 2 + y 2 ) = ( x, y ) √ x 2 + y 2 Notice that ∇ g does not exist at (0 , 0) . Away from that point, the denominator is ‖ ( x, y ) ‖ 2 , and hence ∇ g has norm one. 1.2. (a) The first columns are trivially linearly independent. Since the third column is zero, the dimension is 2. (b) The first three vectors are independent, so the dimension is 3. (c) We can remove the third vector (0 , 1 , 0) from the first two, leaving (2 , 0 , 9) and (3 , 0 , 1). These are not multiples of each other, so the dimension is 3. (d) The first two columns are repeated, so the dimension is 3. 1.3. (a) Linear: 0 + 0 = 0 · 0 = 0 (b) Not linear: 0 · f (0 , 0 , 0) = 1 , but f (0 · 0 , 0 · 0 , 0 · 0) = 1 (c) Not linear: 0 · f (0 , 0 , 0) = (0 , 0) , but f (0 · 0 , 0 · 0 , 0 · 0) = (1 , 0) (d) Linear: both expressions just scale x (e) Linear: Can be written f ( x, y ) =   2 3 1 0 0 0   ( x y ) 1.4. Take x, y ∈ U 1 ∩ U 2 and c ∈ R . In particular, for i ∈ { 1 , 2 } , this means x, y ∈ U i (individually). Since the U i ’s are subspaces, we know x + y, cx ∈ U i for i ∈ { 1 , 2 } . Since they are in both sets individually, we know x + y, cx ∈ U 1 ∩ U 2 , showing that U 1 ∩ U 2 is a subspace. Suppose U 1 = span { (1 , 0) } and U 2 = span { (0 , 1) } . Then, (1 , 0) , (0 , 1) ∈ U 1 ∪ U 2 but (1 , 1) = (1 , 0) + (0 , 1)  ∈ U 1 ∪ U 2 , so U 1 ∪ U 2 cannot be a subspace. 1.5. Suppose f (  x ) = ‖ A x −  a ‖ 2 2 + ‖ B x −  b ‖ 2 2 . Then, at a minimum  0 = ∇ f (  x ) = ∇ [ ‖ A x −  a ‖ 2 2 + ‖ B x −  b ‖ 2 2 ] = ∇ [  x  A  A x − 2  x  A   a + ‖  a ‖ 2 2 +  x  B  B x − 2  x  B   b + ‖  b ‖ 2 2 ] = 2 A  A x − 2 A   a + 2 B  B x − 2 B   b. After dividing by 2 and reorganizing the sides, this shows ( A  A + B  B )  x = A   a + B   b. 1.6. From elementary calculus, f, g ∈ C 1 ( R ) , c ∈ R = ⇒ f + g, cf ∈ C 1 ( R ) . These identities show C 1 ( R ) is a subspace of the space of functions. C 1 ( R ) contains the monomials m k ( x ) ≡ x k . There are infinitely many of these monomials, and each m k cannot be written as a linear combination of m 0 , . . . , m k − 1 with coefficients in R , showing C 1 ( R ) has dimension ∞ . 1.7. ( A  A ) ij =  c i ·  c j and ( AA  ) ij =  r i ·  r j . 1.8. Notice that minimizing f (  x ) is equivalent to minimizing f (  x ) 2 so long as f (  x ) ≥ 0 ∀  x ∈ R . After applying this transformation, taking B = 0 and  b =  0 in problem 1.5 shows A  A x = A   b. 1.9. As in the previous problem, we can instead minimize ‖ A x ‖ 2 2 subject to ‖ B x ‖ 2 2 = 1 . The Lagrange multiplier of the squared problem is Λ(  x ; λ ) ≡  x  A  A x − λ (  x  B  B x − 1) . Differentiating with respect to  x shows  0 = 2 A  A x − 2 λB  B x. This provides the optimality condition: A  A x = λB  B x. 9781482251944_SM_Cover.indd 7 28/09/15 12:38 pm

  4. Page 4
    Read page 4 as text

    2  At this critical point, we can write ‖ A x ‖ 2 = √ ‖ A x ‖ 2 2 = √  x  A  A x = √ λ √  x  B  B x by the optimality condition = √ λ ‖ B x ‖ 2 = √ λ by the constraint. 1.10. After squaring the constraint, the Lagrange multiplier of this system is Λ(  x ; λ ) ≡  a ·  x − λ ( ‖  x ‖ 2 2 − 1) . Differentiating with respect to  x shows  0 =  a − 2 λ x. Plugging this into the constraint, 1 = ‖  x ‖ 2 2 = ‖  a ‖ 2 2 4 λ 2 = ⇒ λ = ± 1 2 ‖  a ‖ 2 . Hence,  x =  a 2 λ = ±  a ‖  a ‖ 2 . Since we wish to maximize f (  x ), we should take  x =  a / ‖  a ‖ 2 , with f (  x ) = ‖  a ‖ 2 . 1.11. Differentiating with respect to  x shows:  0 = ∇ R (  x ) = ∇ (  x  A x ‖  x ‖ 2 2 ) = 2 ‖  x ‖ 2 2 A x − 2(  x  A x )  x ‖  x ‖ 4 2 by the quotient rule = ⇒ A x = (  x  A x ‖  x ‖ 2 2 )  x. Hence,  x is an eigenvector of A with eigenvalue  x  A x / ‖  x ‖ 2 2 . 1.12. By definition, AA − 1 = I n × n . Taking transposes of both sides shows ( A − 1 )  A  = I n × n . Right-multiplying by ( A  ) − 1 shows ( A − 1 )  = ( A  ) − 1 , as desired. Similarly, ( AB )( AB ) − 1 = I n × n = ⇒ B ( AB ) − 1 = A − 1 = ⇒ ( AB ) − 1 = B − 1 A − 1 . 1.13. For two functions A ( t ) , B ( t ) : R → R n × n , the identity d dt ( A ( t ) · B ( t )) = dA dt B + A dB dt follows directly from the product rule in one variable. By definition of matrix inverses A − 1 ( t ) · A ( t ) = I n × n for all t . The right side is a constant with respect to t , so 0 = d dt ( A − 1 ( t ) · A ( t )) = d ( A − 1 ) dt A + A − 1 dA dt = ⇒ d ( A − 1 ) dt = − A − 1 dA dt A − 1 . 1.14. By the chain rule, d dt f (  x + t v ) = ∇ f (  x + t v ) · d dt (  x + t v ) =  v · ∇ f (  x + t v ) . Substituting t = 0 provides the desired result. 9781482251944_SM_Cover.indd 8 28/09/15 12:38 pm

  5. Page 5
    Read page 5 as text

     3 1.15. (a) We expand the product: [ B ( B  B ) − 1 B  ] 2 = [ B ( B  B ) − 1 B  ] · [ B ( B  B ) − 1 B  ] = B [( B  B ) − 1 ( B  B )( B  B ) − 1 ] B  by associativity = B ( B  B ) − 1 B  by canceling out inverses. (b) Again expanding the product, ( I n × n − A ) 2 = I n × n − 2 A + A 2 = I n × n − 2 A + A since A is idempotent = I n × n − A (c) The inverse of 1 2 I n × n − A is 2 I n × n − 4 A , since: ( 1 2 I n × n − A ) (2 I n × n − 4 A ) = I n × n − 2 A − 2 A + 4 A 2 = I n × n − 2 A − 2 A + 4 A since A is idempotent = I n × n (d) A x = λ x = ⇒ A 2  x = λA x = λ 2  x. But A 2 = A , so A 2  x = A x = λ x. Thus, λ x = λ 2  x = ⇒ λ = λ 2 = ⇒ λ ∈ { 0 , 1 } . 1.16. The product AB is of size n 2 . Obviously finding AB requires considering each element of AB at least once (if nothing else, to write the result in memory!), already requiring O ( n 2 ) time even if each element of AB is computed in O (1) time. The algorithms in the figure take O ( n 3 ) time to run due to the nested loops. Hence, there is room for improvement, and indeed Strassen’s algorithm and several others achieve faster than O ( n 3 ) asymptotic runtime, at least for large n . 1.17. Define  (  x ) ≡ − ln p (  x ) . Since ln is monotonic, any local maximum of p (  x ) is also a local maximum of  (  x ) . Hence,  x ∗ is a critical point of  (  x ) , implying ∇  (  x ∗ ) =  0 . Let H be the Hessian of  at  x ∗ . Then, near  x ∗ we can approximate: − ln p (  x ) =  (  x ) ≈  (  x ∗ ) + 1 2 (  x −  x ∗ )  H (  x −  x ∗ ) = − ln p (  x ∗ ) + 1 2 (  x −  x ∗ )  H (  x −  x ∗ ) . The first derivative term of the expansion vanishes since ∇  (  x ∗ ) =  0 . Exponentiating both sides shows p (  x ) ≈ const . · e − 1 2 (  x −  x ∗ )  ( − H )(  x −  x ∗ ) . Hence, a reasonable Gaussian approximation of p (  x ) near  x ∗ takes Σ = − H − 1 and  μ =  x ∗ . 2.1. Depending on the processor, fixed-point arithmetic can be faster than floating-point since it can be carried out on the ALU with integer-type operations without the need for dealing with an exponent. Fixed-point arithmetic also can be applicable when the scale of numbers under consideration is known ahead of time, e.g., in financial software. Floating-point representations are more accurate, especially when values care on many scales. 2.2. (a) (answers may vary) Rounding error can come from multiplication and division to find the value n from the other variables. Discretization error can come from the representations of values from the sensors. Modeling error can result from inaccuracies of the Ideal Gas Law and/or failure to account for secondary factors like sensor noise or pollutants. Input error can result from using an inaccurate value of the constant R . (b) From the ideal gas law, we can write n = P V RT . Then, if we measure ¯ P and ¯ T rather than the ground-truth values, we can write the forward error as: ∣ ∣ ∣ ∣ ¯ P V R ¯ T − P V RT ∣ ∣ ∣ ∣ = V R ∣ ∣ ∣ ∣ ¯ P ¯ T − P T ∣ ∣ ∣ ∣ 9781482251944_SM_Cover.indd 9 28/09/15 12:38 pm

  6. Page 6
    Read page 6 as text

    4  = V R ∣ ∣ ∣ ∣ P + δ P T + δ T − P T ∣ ∣ ∣ ∣ for | δ P | ≤ ε P , | δ T | ≤ ε T = V R ∣ ∣ ∣ ∣ ( P + δ P ) T − P ( T + δ T ) T ( T + δ T ) ∣ ∣ ∣ ∣ = V RT ∣ ∣ ∣ ∣ P T + δ P T − P T − δ T P T + δ T ∣ ∣ ∣ ∣ = n ∣ ∣ ∣ ∣ δ P T − P δ T P ( T + δ T ) ∣ ∣ ∣ ∣ Hence, the relative forward error can be bounded as follows: ∣ ∣ ∣ ∣ δ P T − P δ T P ( T + δ T ) ∣ ∣ ∣ ∣ ≤ T ε P + P ε T P ( T − ε T ) (c) In this case, n = P V RT = (100 Pa)(0 . 5 m 3 ) (8 . 31 J · mol − 1 · K − 1 )(300 K) = 0 . 0201 mol With the given measurement bounds, the largest possible value is (101 Pa)(0 . 5 m 3 ) (8 . 31 J · mol − 1 · K − 1 )(299 . 5 K) = 0 . 0203 mol = n + 0 . 000234 mol = n + 1 . 17% . The smallest possible value is (99 Pa)(0 . 5 m 3 ) (8 . 31 J · mol − 1 · K − 1 )(300 . 5 K) = 0 . 0198 mol = n − 0 . 000234 mol = n − 1 . 17% . Hence, the absolute error is bounded by 0 . 0198 mol and the relative error is bounded by 1 . 17% . (d) At the range indicated by the problem, it is relatively well-conditioned. When the scale of ε T is commen- surate with that of T , the problem becomes ill-conditioned. 2.3. We can understand the relative error as the fraction κ rel = | ∆ y | / | y | | ∆ x | / | x | = ∣ ∣ ∣ ∣ x ∆ y y ∆ x ∣ ∣ ∣ ∣ , where y + ∆ y = f ( x + ∆ x ) and y = f ( x ) . By Taylor’s theorem, f ( x + ∆ x ) = y + f ′ ( x )∆ x + O (∆ x 2 ) . Hence, ∆ y = f ′ ( x )∆ x + O (∆ x 2 ) , so for small ∆ x , κ rel ≈ ∣ ∣ ∣ ∣ x · f ′ ( x )∆ x f ( x ) · ∆ x ∣ ∣ ∣ ∣ = ∣ ∣ ∣ ∣ xf ′ ( x ) f ( x ) ∣ ∣ ∣ ∣ . The absolute condition number of this problem is: ∣ ∣ ∣ ∣ ∆ y ∆ x ∣ ∣ ∣ ∣ ≈ | f ′ ( x ) | . The function f ( x ) = ln x has a large relative condition number near x = 1, since κ rel = 1 / ln x , which blows up near x = 1. Contrastingly, the function f ( x ) = x has relative condition number 1 for all x . 2.4. Since minima are roots of f ′ , we can use the conditioning for root-finding, but with an extra derivative: (a) | x est − x ∗ | (b) | f ′ ( x est ) − f ′ ( x ∗ ) | ≈ δx | f ′′ ( x ∗ ) | (c) 1 / | f ′′ ( x ∗ ) | 2.5. (a) The range is ( −∞ , 0] since lim t → 0 log t = −∞ and log 1 = 0. (b) If the x k is very negative, then e x k is exponentially close to zero. This near-zero value may not be repre- sentable, and regardless a single slightly larger value will dominate the sum. 9781482251944_SM_Cover.indd 10 28/09/15 12:38 pm

  7. Page 7
    Read page 7 as text

     5 (c) We simplify directly:  ( x 1 , . . . , x n ) = ln [ ∑ k e x k ] by definition = ln [ ∑ k e x k − a + a ] = ln [ e a ∑ k e x k − a ] = ln e a + ln [ ∑ k e x k − a ] = a + ln [ ∑ k e x k − a ] Suppose we take a = min k x k . Then, rather than adding together tiny values we have moved the scale to be around e 0 = 1 . (Other heuristics for choosing a are possible) 2.6. There are rendering artifacts because the two surfaces overlap and hence have the same depth values; rounding during depth computation can make one surface appear on top of the other. Possible resolutions include slightly offsetting one surface, adding a tie-breaking rule when depths are within some tolerance of each other, or merging the geometry before rendering to avoid overlap altogether. 2.7. (a) Recall that floating point arithmetic changes spacing as the order of magnitude of the value changes. Thus, it makes sense to have multiplicative error that is relative to the scale of x and y . (b) (adapted from course notes by D. Bindel, Cornell CS) The recurrence for the ground-truth sum is simply s k = s k − 1 + x k y k . Error terms for the addition and multiplication steps show ˆ s k = (ˆ s k − 1 + x k y k (1 + ε × k ))(1 + ε + k ) . Subtracting the two shows: ˆ s k − s k = [(ˆ s k − 1 + x k y k (1 + ε × k ))(1 + ε + k )] − [ s k − 1 + x k y k ] by the recurrences above = [ˆ s k − 1 + ε + k ˆ s k − 1 + x k y k (1 + ε × k ) + x k y k ε + k (1 + ε × k )] − [ s k − 1 + x k y k ] = [ˆ s k − 1 − s k − 1 ] + ε + k ˆ s k − 1 + x k y k ( ε × k + ε + k + ε + k ε × k ) = [ˆ s k − 1 − s k − 1 ](1 + ε + k ) + ε + k s k − 1 + x k y k ( ε × k + ε + k + ε + k ε × k ) = [ˆ s k − 1 − s k − 1 ](1 + ε + k ) + ε + k s k + x k y k ( ε × k + ε + k ε × k ) since s k = s k − 1 + x k y k = [ˆ s k − 1 − s k − 1 ](1 + ε + k ) + ε + k s k + x k y k ε × k + x k y k ε + k ε × k We can expand this inductively: ˆ s 0 − s 0 = 0 ˆ s 1 − s 1 = [ˆ s 0 − s 0 ](1 + ε + 1 ) + ε + 1 x 1 y 1 + x 1 y 1 ε × 1 + x 1 y 1 ε + 1 ε × 1 = x 1 y 1 ( ε + 1 + ε × 1 ) + x 1 y 1 ε + 1 ε × 1 ˆ s 2 − s 2 = [ˆ s 1 − s 1 ](1 + ε + 2 ) + ε + 2 ( x 1 y 1 + x 2 y 2 ) + x 2 y 2 ε × 2 + x 2 y 2 ε + 2 ε × 2 = [ x 1 y 1 ( ε + 1 + ε × 1 ) + x 1 y 1 ε + 1 ε × 1 ](1 + ε + 2 ) + ε + 2 ( x 1 y 1 + x 2 y 2 ) + x 2 y 2 ε × 2 + x 2 y 2 ε + 2 ε × 2 = x 1 y 1 ( ε + 1 + ε × 1 + ε + 1 ε + 2 + ε × 1 ε + 2 + ε + 2 ) + x 2 y 2 ( ε + 2 + ε × 2 ) + [ x 1 y 1 ε + 1 ε × 1 + x 2 y 2 ε + 2 ε × 2 ] + O ( ε 3 max ) . . . Applying induction, this recurrence shows ˆ s k − s k = k ∑ i =1   x i y i   ε × i + k ∑ j = i ε + j     + O ( kε 2 max ) = ⇒ e n ≤ nε max ∑ k | x k || y k | + O ( nε 2 max ) , as desired. 9781482251944_SM_Cover.indd 11 28/09/15 12:38 pm

  8. Page 8
    Read page 8 as text

    6  2.8. For convenience, define d ≡ x − y . We’ll start by simplifying the numerator of relative error and then substitute: (1 + ε x ) x − (1 + ε y ) y = ( x − y ) + ( ε x x − ε y y ) = d + ε x d + ( ε x − ε y ) y = ⇒ (1 + ε − )((1 + ε x ) x − (1 + ε y ) y ) = (1 + ε − )( d + ε x d + ( ε x − ε y ) y ) = (1 + ε − ) d + ε x (1 + ε − ) d + (1 + ε − )( ε x − ε y ) y = ⇒ E = ∣ ∣ ∣ ∣ (1 + ε − )((1 + ε x ) x − (1 + ε y ) y ) − ( x − y ) x − y ∣ ∣ ∣ ∣ = ∣ ∣ ∣ ∣ ε − d + ε x (1 + ε − ) d + (1 + ε − )( ε x − ε y ) y d ∣ ∣ ∣ ∣ = ∣ ∣ ∣ ε − + ε x (1 + ε − ) + (1 + ε − )( ε x − ε y ) y d ∣ ∣ ∣ This can be unbounded as d → 0. 2.9. (a) Implicitly differentiating the relationship 0 = f ( x ( ε )) + εp ( x ( ε )) with respect to ε shows 0 = d dε [ f ( x ( ε )) + εp ( x ( ε ))] = f ′ ( x ( ε )) x ′ ( ε ) + p ( x ( ε )) + εp ′ ( x ( ε )) x ′ ( ε ) by the chain rule. Substituting ε = 0 and using x ∗ = x (0) shows 0 = f ( x ∗ ) x ′ (0) + p ( x ∗ ) = ⇒ x ′ (0) = − p ( x ∗ ) f ( x ∗ ) . (b) We differentiate f ′ ( x ) = d dx ( x − 1) · ( x − 2) · · · · · ( x − 20) = ( x − 2) · · · · · ( x − 20) + ( x − 1) · ( x − 3) · · · · · ( x − 20) + · · · + ( x − 1) · · · · · ( x − 19) by the product rule Substituting x = j shows f ′ ( j ) = ( j − 1) · ( j − 2) · · · · · ( j − ( j − 1)) · ( j − ( j + 1)) · · · · · ( j − 20) For p ( x ) = x 19 , from the previous part we have x ′ ( j ) = − j 19 ( j − 1) · ( j − 2) · · · · · ( j − ( j − 1)) · ( j − ( j + 1)) · · · · · ( j − 20) = − ∏ k  = j j j − k . (c) x ′ (1) ≈ 8 . 2 × 10 − 18 and x ′ (20) ≈ − 4 . 3 × 10 7 ; hence, the root x ∗ = 1 is far more stable. 2.10. (a) The alternative formula can be obtained by scaling the numerator and denominator of the quadratic equation: − b ± √ b 2 − 4 ac 2 a = − b ± √ b 2 − 4 ac 2 a · − b ∓ √ b 2 − 4 ac − b ∓ √ b 2 − 4 ac = b 2 − ( b 2 − 4 ac ) − 2 ab ∓ 2 a √ b 2 − 4 ac = 4 ac − 2 ab ∓ 2 a √ b 2 − 4 ac = − 2 c b ± √ b 2 − 4 ac 9781482251944_SM_Cover.indd 12 28/09/15 12:38 pm

  9. Page 9
    Read page 9 as text

     7 (b) When b ≤ 0, take x 1 = − b + √ b 2 − 4 ac 2 a , x 2 = c ax 1 , and otherwise take x 2 = − b − √ b 2 − 4 ac 2 a , x 2 = c ax 2 . This way, there never can be cancellation because we always move b farther from the origin in the numer- ator. 2.11. The bounds are worked out below: [ x ] + [ y ] = [ x + y, x + y ] [ x ] − [ y ] = [ x − y, x − y ] [ x ] × [ y ] =                                value sign( x ) sign( x ) sign( y ) sign( y ) [ xy , x y ] + + + + [ xy , x y ] + + − + [ y x, yx ] + + − − [ x y, xy ] − + + + [ min( x y, y x ) , max( xy, xy ) ] − + − + [ xy, xy ] − + − − [ x y, xy ] − − + + [ x y, x y ] − − − + [ x y, xy ] − − − − [ x ] ÷ [ y ] = [ x ] × [ 1 y , 1 y ] [ x ] 1 / 2 = [ x 1 / 2 , x 1 / 2 ] In finite-precision arithmetic, always round down the lower bounds and round up the upper bounds. 2.12. (a) Perturbing any of three collinear points slightly makes them not collinear. Furthermore, points may appear collinear if you zoom out far enough but appear less so as you zoom in. (b) ε ε ε ε ε ε  p  q 9781482251944_SM_Cover.indd 13 28/09/15 12:38 pm

  10. Page 10
    Read page 10 as text

    8  (c) ε ε ε ε  p  q (d) Obvious from drawings above; ε -collinear points form the intersection of four half-planes, two of which come from the ε -clockwise condition and two of which come from the ε -counterclockwise condition. (e) No. See § 3.1 of [55] for an example. 3.1. No; LU may not be possible for matrices requiring pivoting. 3.2. The steps of Gaussian elimination are below: ( 2 4 2 3 5 4 ) ∼ ( 1 2 1 3 5 4 ) , with elimination matrix ( 1 / 2 0 0 1 ) ∼ ( 1 2 1 0 1 − 1 ) , with elimination matrix ( 1 0 3 − 1 ) ∼ ( 1 0 3 0 1 − 1 ) , with elimination matrix ( 1 − 2 0 1 ) So, x = 3 and y = − 1 . From the steps above, we know U = ( 1 2 0 1 ) , and L = ( 1 / 2 0 0 1 ) − 1 ( 1 0 3 − 1 ) − 1 = ( 2 0 0 1 ) ( 1 0 3 − 1 ) = ( 2 0 3 − 1 ) . 3.3. Computed using Gaussian elimination: L =   1 0 0 3 1 0 6 11 1   U =   1 2 7 0 − 1 − 22 0 0 204   3.4. Where it states “optionally insert pivoting code here,” find row r with largest value in column p ; then swap row r and row p of both A and  b. 3.5. No. Full pivoting can be preferable numerically but technically does not make a difference. The only way partial pivoting would fail is if there is an all-zero column, which would indicate that A is not invertible. 3.6. Write A = A 1 + A 2 i ,  b =  b 1 +  b 2 i , and  x =  x 1 +  x 2 i. Then, A x =  b = ⇒ ( A 1 + A 2 i )(  x 1 +  x 2 i ) =  b 1 +  b 2 i = ⇒ ( A 1  x 1 − A 2  x 2 ) + ( A 2  x 1 + A 1  x 2 ) i =  b 1 +  b 2 i. So, we can solve the block system ( A 1 − A 2 A 2 A 1 ) (  x 1  x 2 ) = (  b 1  b 2 ) . 3.7. Carrying out Gaussian elimination is the same as pre-multiplying by the inverse of the leftmost n × n block. Hence, the output is A − 1 ( A | I n × n ) = ( A − 1 A | A − 1 ) = ( I n × n | A − 1 ) . 9781482251944_SM_Cover.indd 14 28/09/15 12:38 pm

  11. Page 11
    Read page 11 as text

     9 3.8. The determinant of a lower triangular matrix is the product of its diagonal elements. If L has a zero diagonal element, det( L ) = 0, contradicting that L is invertible. We used this in § 3.5.3 when we assumed L can have ones along the diagonal. 3.9. Lower triangular matrices are solved using forward substitution. Forward substitution matrices are lower- triangular, so the inverse of a lower triangular matrix is lower-triangular. By Proposition 3.1, the inverse is thus lower-triangular. 3.10. Notice a 11 =  11 u 11 . So if a 11 = 0 then one of  11 or u 11 must equal zero. This trivially implies that L or U isn’t full rank, contradicting A being full rank. 3.11. From linear algebra, det( A ) = det( L ) · det( U ) , and determinant of a triangular matrix is the product of its diagonal elements. 3.12. (a) The initial forward substitution step yields an upper triangular matrix U . U is obtained by alternating be- tween pivoting (a permutation matrix) and forward substitution (a product of lower-triangular elimination matrices). (b) P i is the matrix of partial pivoting for row i ; partial pivoting only considers unprocessed rows j ≥ i . Similarly, each forward substitution matrix takes the form indicated because it replicates a multiple of row i to a row k > i ; see § 3.3.3. (c) We simplify as follows: P jk ( I n × n + c e k  e  i ) = P jk + c ( P jk  e k )  e  i by distributing the product = P jk + c e j  e  i since P jk permutes rows j and k = P jk + c e j ( P jk  e i )  since P jk does not affect row i = P jk + c e j  e  i P  jk since ( AB )  = B  A  = P jk + c e j  e  i P jk since the inverse of a row swap is the same row swap = ( I n × n + c e j  e  i ) P jk , as needed. (d) Follows from induction. Parts (b) and (c) together show that P i L j = L ′ P i for some lower-triangular L ′ when i > j . (e) From (a), U = L n − 1 P n − 1 · · · L 1 P 1 A . Applying (d) shows U = LP A for some lower-triangular matrix L and some permutation matrix P . Hence, A = P  L − 1 U, and the result follows from the facts that the transpose of a permutation is a permutation and the inverse of a lower-triangular matrix is lower-triangular. (f) Suppose A = P LU. Then, A − 1  b = U − 1 L − 1 P   b. So, apply P  to  b , then L − 1 via forward substitution, then U − 1 via back substitution. 3.13. (a) We can work through the product in the reverse direction: ( I 0 CA − 1 I ) ( A 0 0 D − CA − 1 B ) ( I A − 1 B 0 I ) = ( I 0 CA − 1 I ) ( A B 0 D − CA − 1 B ) = ( A B C D ) = M , as desired. (b) In our factorization from (a), the left-hand and right-hand matrices are lower- and upper-triangular, respectively. Hence, we only have to factor the middle matrix as a product LU to get our factorization (the product of triangular matrices is triangular). We expand: ( A 0 0 D − CA − 1 B ) = ( L 1 U 1 0 0 L 2 U 2 ) , as constructed = ( L 1 0 0 L 2 ) ( U 1 0 0 U 2 ) Hence, from the factorization in (a), we can write: M = [( I 0 CA − 1 I ) ( L 1 0 0 L 2 )]  L [( U 1 0 0 U 2 ) ( I A − 1 B 0 I )]  U 9781482251944_SM_Cover.indd 15 28/09/15 12:38 pm

  12. Page 12
    Read page 12 as text

    10  (c) The recursive algorithm could proceed as follows (runtime denoted f ( n ) for factoring an n × n matrix): • Recursively factor A = L 1 U 1 (runtime f ( n / 2 )). • Construct A − 1 B = U − 1 1 L − 1 1 B using forward and then backward substitution (runtime O ( n 3 )). • Construct R ≡ D − CA − 1 B via matrix multiplication/subtraction (runtime O ( n 3 )). • Factor R = L 2 U 2 (runtime f ( n / 2 )). • Return the factorization from (b) (runtime O ( n 3 ) for matrix multiplication and construction of CA − 1 ). Hence, the runtime satisfies f ( n ) = 2 f ( n / 2 ) + O ( n 3 ), which is slower than the O ( n 3 ) method presented in the chapter. As an aside, however, a few recursive steps using this formula may help break the LU factorization problem up for parallel algorithmic processing. 3.14. It suffices to show that forward substitution can be carried out without pivoting, each time leaving a lower right rectangle that is diagonally dominant. This ensures that pivots are always nonzero, and subsequently back substitution can be carried out as usual. Suppose we write A as A = ( c  v   w B ) . After Gaussian elimination, we modify to ¯ A = ( 1 1 c  v   0 B −  w v  c ) . So, we must check that ¯ B ≡ B −  w v  c is diagonally dominant. We do so for column i as follows: ∑ j  = i | ¯ b ji | = ∑ j  = i | b ji − w j v i /c | ≤ ∑ j  = i | b ji | + | v i | | c | ∑ j  = i | w j | by the triangle inequality ≤ ( | b ii | − | v i | ) + | v i | | c | ∑ j  = i | w j | by diagonal dominance of A in column i ≤ ( | b ii | − | v i | ) + | v i | | c | ( | c | − | w i | ) by diagonal dominance of A in column 1 = | b ii | − | v i | | c | | w i | ≤ | b ii − w i v i /c | by the triangle inequality = | ¯ b ii | , as needed. From a numerical standpoint, pivoting may help improve conditioning of arithmetic operations, as dividing by c can still cause issues. 3.15. Suppose A = LU = L ′ U ′ , where L and L ′ have ones along their diagonals. Then, ( L ′ ) − 1 L = U ′ U − 1 . The left-hand side of this equality is lower triangular and the right-hand side is upper triangular. Hence, D ≡ ( L ′ ) − 1 L = U ′ U − 1 is a diagonal matrix. Writing L = L ′ D shows D = I n × n since L, L ′ have ones along their diagonals (and D scales L ′ column-by-column). In other words, L = L ′ , as needed. 4.1. Many possible answers. One example is the matrix whose diagonal and superdiagonal both equal 1, e.g., for size 6 × 6: A =         1 1 0 0 0 0 0 1 1 0 0 0 0 0 1 1 0 0 0 0 0 1 1 0 0 0 0 0 1 1 0 0 0 0 0 1         . 9781482251944_SM_Cover.indd 16 28/09/15 12:38 pm

  13. Page 13
    Read page 13 as text

     11 This matrix has an upper-triangular inverse that alternates between ± 1: A − 1 =         1 − 1 1 − 1 1 − 1 0 1 − 1 1 − 1 1 0 0 1 − 1 1 − 1 0 0 0 1 − 1 1 0 0 0 0 1 − 1 0 0 0 0 0 1         4.2. D is the lower-right submatrix of ECE  ; since ECE  is symmetric (by symmetry of C ), D must be symmetric as well. Positive-definiteness of D follows from positive-definiteness of C and invertibility of E :  w  D  w = ( 0  w )  ECE  ( 0  w ) > 0 . 4.3. (a) The gradient of f contains the partial derivatives with respect to each of n × m possible perturbations, which can be organized naturally in an n × m matrix. (b) (i)  x y  (ii) 2 A x x  (iii) − 2 W (  x − A y )  y  (c) (i) Follows from linearity of d / dt (ii) Trivial since ∂ is taken elementwise; the transpose simply re-indexes the elements of X (iii) We work this one out more explicitly: ( XY ) ′ ij = ∑ k ( X ik Y kj ) ′ by definition of matrix multiplication = ∑ k ( X ′ ik Y kj + X ik Y ′ kj ) by the (scalar) product rule = ( ∂X ) Y + X ( ∂Y ) by definition of matrix multiplication (iv) Since XX − 1 = I n × n , form the previous part we know 0 = ( ∂X ) X − 1 + X ( ∂X − 1 ) . Isolating ∂X − 1 provides the desired result. 4.4. Define  x ≡  x sum ,  y ≡  y sum . We know A x + p  b =  y. More explicitly, a 11 x 1 + a 12 x 2 + pb 1 = y 1 a 21 x 1 + a 22 x 2 + pb 2 = y 2 Similarly, AX +  b x  = C implies: a 11 X 11 + a 12 X 21 + b 1 x 1 = C 11 a 11 X 12 + a 12 X 22 + b 1 x 2 = C 12 a 21 X 11 + a 22 X 21 + b 2 x 1 = C 21 a 21 X 12 + a 22 X 22 + b 2 x 2 = C 22 Hence, we can write:         x 1 x 2 0 0 p 0 0 0 x 1 x 2 0 p X 11 X 21 0 0 x 1 0 X 12 X 22 0 0 x 2 0 0 0 X 11 X 21 0 x 1 0 0 X 12 X 22 0 x 2          M         a 11 a 12 a 21 a 22 b 1 b 2           u =         y 1 y 2 C 11 C 12 C 21 C 22           d . The systems for x and y components decouple, as might be expected. 9781482251944_SM_Cover.indd 17 28/09/15 12:38 pm

  14. Page 14
    Read page 14 as text

    12  4.5. Rather than using Lagrange multipliers, we can simply assume that  t ( v ) is a constant in E D if v ∈ B. Suppose u  ∈ B . Then, 0 = ∂ ∂ t ( u ) ∑ ( v,w ) ∈ E ‖  t ( v ) −  t ( w ) ‖ 2 2 at a critical point = ∑ v ∈ N ( u ) (  t ( u ) −  t ( v )) since all other terms are zero = | N ( u ) |  t ( u ) − ∑ v ∈ N ( u )  t ( v ) since the first term is a constant = ⇒  t ( u ) = 1 | N ( u ) | ∑ v ∈ N ( u )  t ( v ), as desired. 4.6. (a) Suppose A = LU , and define D = U L − . Note LDL  = L ( U L − ) L  = LU = A , so we only need to show that D is diagonal. Since L is lower-triangular, L  and hence L − are both upper-triangular; this shows D is the product of upper-triangular matrices and thus is upper triangular. Since A is symmetric, LDL  = A = A  = LD  L  = ⇒ D = D  . We have shown that D is a symmetric upper-triangular matrix, which implies that it is diagonal. (b) Suppose we re-define E without using square roots as E ≡ ( 1  0   r I ( n − 1) × ( n − 1) ) . Then, EC = ( c 11  v  0 D ) = ⇒ ECE  = ( c 11  0   0 D ) . We can apply the same pattern of symmetric pre- and post-multiplication, but now get a diagonal matrix down the center rather than the identity. This algorithm did not need square roots and hence can work even when the matrix is symmetric but not positive definite. 4.7. We equivalently minimize ‖  x ‖ 2 2 / 2 rather than ‖  x ‖ 2 . Then, the Lagrange multiplier system becomes Λ = 1 2 ‖  x ‖ 2 2 +  λ  (  b − A x ) = ⇒ 0 = ∇  x Λ =  x − A   λ = ⇒  x = A   λ. Returning to the constraint,  b = A x = AA   λ = ⇒  λ = ( AA  ) − 1  b. Substituting this into the expression above shows,  x = A   λ = A  ( AA  ) − 1  b. Recall that the Tikhonov system minimizes ‖ A x −  b ‖ 2 2 + α ‖  x ‖ 2 2 for some α > 0. We can decompose any  x ∈ R n as  x = A   y +  x 0 , where A x 0 =  0. In this notation, ‖  x ‖ 2 2 =  x   x = ( A   y +  x 0 )  ( A   y +  x 0 ) by our definition =  y  AA   y +  y  A x 0 +  x  0 A   y +  x  0  x 0 by expanding the square =  y  AA   y +  x  0  x 0 since A x 0 =  0 . Similarly, ‖ A x ‖ 2 2 = ‖ A ( A   y +  x 0 ) ‖ 2 2 by our definition = ‖ AA   y ‖ 2 2 since A x 0 =  0  b  A x =  b  AA   y From our substitution and subsequent formulas, an equivalent way to carry out Tikhonov-regularized least- squares is to solve the following optimization: min  x 0 , y α y  AA   y + α x  0  x 0 + ‖ AA   y ‖ 2 2 − 2  b  AA   y 9781482251944_SM_Cover.indd 18 28/09/15 12:38 pm

  15. Page 15
    Read page 15 as text

     13 s.t. A   x 0 =  0 . The problem decouples over  x 0 and  y . The optimal  x 0 is  0. The optimality condition for  y is  0 = αAA   y + ( AA  ) 2  y − AA   b = ⇒  y = (( AA  ) + αI n × n ) − 1  b = ⇒  x = A   y +  x 0 = A  (( AA  ) + αI n × n ) − 1  b +  0 from our work above → A  ( AA  ) − 1  b as α → 0, as desired. The final limit is permissible because AA  is invertible; any argument involving “( A  A ) − 1 ” cannot be correct because A  A is not necessarily invertible. 4.8. The algorithm takes a forward pass and then a backward pass. The forward pass uses these formulas for increasing i : ¯ w i = { w i / v i when i = 1 w i / ( v i − u i ¯ w i − 1 ) when i > 1 ¯ b i = { b i / v i when i = 1 ( b i − u 2 ¯ b i − 1 ) / ( v i − u i ¯ w i − 1 ) when i > 1 . The backward pass becomes: x n = ¯ b n x i = ¯ x i − ¯ w i x i +1 . 4.9. The Lagrange multiplier for this problem is Λ = 1 2  x  A  A x +  λ  (  c − B x ) . Differentiating with respect to  x shows  0 = A  A x − B   λ. Hence, we can solve the following linear system: ( A  A − B  B 0 ) (  x  λ ) = (  0  c ) . 4.10. (a) Take  x ∈ R n = ⇒  x  A x =  x  LL   x = ‖ L   x ‖ 2 2 ≥ 0 . (b) Try to run Cholesky factorization; it is positive semidefinite only if the factorization succeeds. 4.11. Yes; R m × n is a finite-dimensional ( mn -dimensional) vector space, and all norms on finite-dimensional vector spaces are equivalent. 4.12. (a) max ‖  x ‖ =1 ‖ ( A + B )  x ‖ ≤ max ‖  x ‖ =1 ( ‖ A x ‖ + ‖ B x ‖ ) by the triangle inequality for vector norms ≤ ( max ‖  x ‖ =1 ‖ A x ‖ ) + ( max ‖  y ‖ =1 ‖ B y ‖ ) since this is a looser optimization = ‖ A ‖ + ‖ B ‖ (b) ‖ A v ‖ = ‖  v ‖ ∥ ∥ ∥ ∥ A  v ‖  v ‖ ∥ ∥ ∥ ∥ ≤ ‖  v ‖ max ‖  x ‖ =1 ‖ A x ‖ since  v / ‖  v ‖ has unit length = ‖ A ‖‖  v ‖ ‖ AB ‖ = max ‖  x ‖ =1 ‖ AB x ‖ ≤ max ‖  x ‖ =1 ( ‖ A ‖‖ B x ‖ ) by the last part = ‖ A ‖ max ‖  x ‖ =1 ‖ B x ‖ = ‖ A ‖‖ B ‖ 9781482251944_SM_Cover.indd 19 28/09/15 12:38 pm

  16. Page 16
    Read page 16 as text

    14  (c) Take  v to satisfy A v = λ v . Then, | λ | k ‖  v ‖ = ‖ A k  v ‖ ≤ ‖ A k ‖‖  v ‖ by a previous part. Dividing by ‖  v ‖ and taking the k -th root provides the desired result. (d) Define the column sum c j = ∑ i | a ij | . We simplify: ‖ A ‖ = max ‖  x ‖ 1 =1 ‖ A x ‖ 1 = max ∑ i | x i | =1 ∑ i ∣ ∣ ∣ ∣ ∣ ∣ ∑ j a ij x j ∣ ∣ ∣ ∣ ∣ ∣ by definition ≤ max ∑ i | x i | =1 ∑ i ∑ j | a ij || x j | by the triangle inequality = max ∑ i | x i | =1 ∑ j | x j | c j Take k ≡ arg max j c j , and without loss of generality since this expression includes absolute values we can assume x j ≥ 0 for all j . Suppose x j > 0 for some j  = k . Then, define  y =  x − x j  e j + x j  e k . We have ∑ j | y j | c j ≥ ∑ j | x j | c j . This shows  x =  e k , and thus ‖ A ‖ ≤ c k . A e k = c k , so this upper bound is achieved. (e) This is a relatively standard theorem; see a matrix analysis text. 4.13. (a) (i) E ( i ) (  z ) ≡ ‖  y −  z ‖ 2 2 (ii) E ( ii ) (  z ) ≡ ∑ i ( z i +1 − z i ) 2 (b) Minimize E ( i ) (  z ) + αE ( ii ) (  z ); all terms are quadratic in elements of  z , so when we take the gradient we will get a linear system (c) Sparse, positive (semi-)definite 4.14. (a) ( XX  + αI )  a = X y (b) Decompose  a =  a X +  a ⊥ where  a X ∈ span {  x ( i ) } and  a ⊥ ·  x ( i ) = 0 for all i . Then, X y = ( XX  + αI )  a = ( XX  + αI )(  a X +  a ⊥ ) = ( XX  + αI )  a X + α a ⊥ by definition of  a ⊥ = ⇒ X (  y − X   a X ) − α a X = α a ⊥ Notice that the left-hand side is in span {  x ( i ) } since the first term is a product with X and by definition of  a X . But, since α  = 0, by definition of  a ⊥ we must have  a ⊥ =  0 since otherwise it would have a component in span {  x ( i ) } . (c) Notice that  a = X c . Thus, we can write ( XX  + αI ) X c = X y . Alternatively, we can factor X ( X  X + αI )  c = X y . Pre-multiplying by X  shows ( X  X )( X  X + αI )  c = ( X  X )  y , and by invertibility of X  X we know ( X  X + αI )  c =  y. (d) Take the columns of X φ to be φ (  x ( i ) ) , and take  c to be the solution of ( X  φ X φ + αI )  c =  y ; notice that the elements of X  φ X φ are dot products φ (  x ( i ) ) · φ (  x ( j ) ) = K (  x ( i ) ,  x ( j ) ), which can be computed without direct knowledge of φ . Then, we can write: f φ (  x ) =  a · φ (  x ) = ( X φ  c ) · φ (  x ) = ∑ i c i φ (  x ( i ) ) · φ ( x ) = ∑ i c i K (  x ( i ) ,  x ) (e) The integral is a generalized inner product between functions of s and t where the sum has been replaced with an integral over x with weights e − πx 2 . The algorithm we developed only needs values K (  x ( i ) ,  x ( j ) ) to compute the c i ’s, and then evaluates f φ as a sum over K (  x ( i ) ,  x ) . Hence, we can carry out the algorithm without ever doing the inner product integral explicitly. 9781482251944_SM_Cover.indd 20 28/09/15 12:38 pm

40 more pages in the full document

Unlock it once and it stays in your library, ready to chat with or turn into flashcards.

Unlock the full document

More Electrical Engineering & Electronics study materials

See all

Study this document with CramX

Once it is in your library, every tool can work from it.