1. 1 BASIC PC BASIC BASIC BASIC Fortran WS PC (1.3) 1 + x 1 x = x = (1.1) 1 + x = (1.2) 1 + x 1 = (1.

Section Title Pages Id 1 3 7239 2 4 7239 3 10 7239 4 8 7244 5 13 7276 6 14 7338 7 8 7338 8 7 7445 9 11 7580 10 10 7590 11 8 7580 12 6 7395 13 z 11 7746 14 13 7753 15 7 7859 16 8 7942 17 8 Id URL http://km.int.oyo.co.jp/showdocumentdetailspage.jsp?documentid= 1

1. 1 BASIC PC BASIC BASIC BASIC Fortran WS PC 1 6 3 (1.3) 1 + x 1 x = 1.23456 10 3 6 1 + x = 1.00123 (1.1) 1 + x = 1.00061 (1.2) 1 + x 1 = 0.00061 (1.3) (1.3) 4 x 1 + x 1 = 2 1 + x + 1 1.23456 10 3 = = 6.17092 10 4 2.00061 (1.1) 1 + x = 1.00123456

1. 2 cos θ > 0 1 cos θ 1 cos θ = sin2 θ θ 0 = 35 θ = 35 30 1 + cos θ x 2 6.28318x + 0.123456 = 0 x 1, x 2 6 = 39.4784 0.493824 = 38.9846 x 1, x 2 = 6.28318 ± 38.9846 2 6.28318 ± 6.24376 = = 6.26347, 0.01971 2 x 1 x 2 s 2 2 = 1 N x 2 i x 2 (1.6) N x 2 = 0.123456 x 1 = 0.0197105 x 2 y 2 x = 37507 y = 37496 x 2 y 2 = 1.40678 10 9 1.40595 10 9 = 8.3 10 5 4 sin θ = 0.580703 sin θ 0 = 0.573576 2 θ θ 0 (1.4) N x i (1 i N) x s 2 x = 1 N N x i (1.5) (x y)(x + y) = 11 75003 = 825033 N s 2 = 119479 (345.613) 2 x y 3 = 119479 119448 = 31 x, y x 2, y 2 10 x s 2 4 x i x sin θ sin θ 0 x = 1 N (x i x) 2 + x (1.7) N i=1 sin θ sin θ 0 = 2 cos θ + θ 0 sin θ θ 0 s 2 = 1 N (x (1.4) i x) 2 ( x x) 2 (1.8) N 2 2 i=1 i=1 i=1 x i i x i i x i i x i 1 348.200 5 350.441 9 352.683 2 338.541 6 335.953 10 346.906 3 355.271 7 343.024 11 344.318 4 342.076 8 349.148 12 340.783 x = 4147.35 = 345.613 12 1 x 2 i = 119479

1. 3 (1.6) x = 0 sin 8 ϕ x = x x x γ µgal 9 335 ϕ 355 ϕ 0 x = 345 γ 0 γ γ 0 = 5.163005(sin 2 ϕ sin 2 ϕ 0 ) i x i x i x i x i x i x [(1 + 0.0088885 sin 2 ϕ 0 ) 1 3.200 5 5.441 9 7.683 +0.0044442(sin 2 ϕ sin 2 ϕ 0 )] 2 6.459 6 9.047 10 1.906 3 10.271 7 1.976 11 0.682 4 2.924 8 4.148 12 4.217 sin 2 ϕ sin 2 ϕ 0 = sin(ϕ + ϕ 0 ) sin(ϕ ϕ 0 ) ϕ ϕ (1.7) (1.8) 0 = ±10 x x = 1 γ γ 0 (xi x) = 0.612 N 5.16 sin 70 sin 10 = 0.84 [gal] s 2 = 379.474 0.612 2 12 = 31.6228 0.374544 = 31.2483 γ γ 0 µgal 6 ϕ 0 = 35 (1/N) x i sin 2 35 = 0.3289899 x (1/N) x 2 i x2 x γ γ x i x 0 = 5.178103 sin(ϕ + 35 ) sin(ϕ 35 ) x [1 + 0.002257 sin(ϕ + 35 ) sin(ϕ 35 )] e x e x = 1 + x 1! + x2 2! + x3 3! + 13 24 49.2 x x = 10 n a n S n ϕ ( ) 15 7.64717e+2 2.98494e+2 γ = 978.032681(1 + 0.005278970 sin 2 16 4.77948e+2 1.79454e+2 ϕ 17 2.81146e+2 1.01692e+2 +0.000023461 sin 4 ϕ) [gal] 18 1.56192e+2 5.45001e+1 4 µgal 19 8.22064e+1 2.77063e+1 20 4.11032e+1 1.33969e+1

2. 4 a n n S n n e+2 10 2 6 654321 0.1 x = +10 n a n S n 15 7.64717e+2 2.09529e+4 16 4.77948e+2 2.14308e+4 17 2.81146e+2 2.17120e+4 18 1.56192e+2 2.18682e+4 19 8.22064e+1 2.19504e+4 20 4.11032e+1 2.19915e+4 ( 2.20265e+4 0.1 10 654321 654322 2 0 1 8 16 4 ) e 10 e 10 1 0 1/e 10 15 16 1 16 10 a b c, d, e, f, 10 0 1 2 3 4 5 6 7 2 0 1 10 11 100 101 110 111 10 8 9 10 11 12 13 14 15 2 1000 1001 1010 1011 1100 1101 1110 1111 8 10 11 12 13 14 15 16 17 16 8 9 a b c d e f n 2 31 = 2147483648 n 2 31 1 = 2147483647 10 2 8-4 100 4-3 101 5-2 110 6 0 000 0 1 001 1 2 2 010 2 3 3 011 3 1 x x = ±(d 1 β 1 + d 2 β 2 + + d t β t )β e 32 = ±mβ e (2.1)

2. 5 0 d i < β 0 m < 1 1 32 β m e (1bit) (8bits) (23bits) 0/1 e + 127 d 1 d 2 d 23 IBM HITACHI 1 β = 16 16 0 1 UNIX 127 = 2 7 1 x = 0.25 1 8 e β = 2 x = 1 2 2 1 e = 1 d 1 = 1 0 2 7 1 + e 2 8 1 127 e 128 m = (0.1000 0) 2 β = 16 x = 4 16 160 e = 0 d 1 = 4 m = (0.0100 00 0) 2 2 127 = 5.9 10 39 2 128 = 3.4 10 38 ( ) 2 16 2.0 = (1.0) 2 2 1 (0100 0000 ) 2 = (4000 0000) 16 x = 0.1 1.5 = (1.1) 2 (0011 1111 1100 0000 ) 2 β = 2 e = 3 m = (0.11001100110 ) 2 = (3fc0 0000) 16 1.0 = (1.0) 0.1 2 (0011 1111 1000 0000 ) 2 = (3f80 0000) 16 Fortran x=0 1 x=x+0.1... if( x.ne.1 ) goto 1... UNIX x = ±(1.d 1 d 2 d 23 ) 2 2 e 0.1 = (1.10011001100 ) 2 2 4 (0011 1101 1100 1100 ) 2 = (3dcc cccd) 16 16 1 16 β 16 (1bit) (7bit) (24bit) UNIX 0/1 e + 64 d 1 d 2 d 24 IEEE 24 UNIX 1 32 1 UNIX 24 16 x

2. 6 16 c/(a b) a b 16 65 = 5.4 10 79 16 63 = 7.2 10 75 0 (machin epsilon) 16 16 ε M 1 36 ε M = Min{ε : 1 + ε > 1} 36 1 + ε > 1 36 6 ε 6 Fortran x 2 + y 2 eps = 1 x y x 2 y 2 do i=1, 64 eps1 = 1 + eps ( ) ( ) x 2 + y 2 if( eps1.eq.1 ) goto 1 eps = eps/2 enddo ( ) stop 1 eps = eps*2 print eps NaN(not a number) eps x y x2 + y 2 = x 1 + (y/x) 2 eps = 5.96046 10 8 c/a/b 1.387778 10 17

2. 7 a + (b + c) = (a + b) + c A/D y = x + x x 0 y = x 1 1 deltax =... ( 16 )... if( y.ne.x+deltax ) goto 1... 16 24 A/D x 0 1 32 16 8 x eps1=1+eps ( ) if(x.eq.y)

3. 1 3 x = u + v (3.5) u 3 + v 3 + q + 3(u + v)(uv + p) = 0 u 3 + v 3 = q uv = p x 2 + 2bx + c = 0 (3.1) u 3 = 1 2 [ q ± q 2 + 4p 3 ] x = b ± b 2 c (3.2) v 3 = p3 u 3 (3.6) u 3 1 b 2 c D = q 2 + 4p 3 (3.7) x 1 = b b 2 c x 2 = c b > 0 x 1 x 1 = b + b 2 c x 2 = c b < 0 (3.3) D u 3 v 3 x 1 u v 1 2 b (3.5) (3.4) b 2 1 2 u 3 v 3 1 y 3 + a 2 y 2 + a 1 x + a 0 = 0 y = x a 2 3 x 3 + 3px + q = 0 (3.4) u 3 = u 1 e 2πi/3 p = 1 3 a 1 [ 2 q = a 2 3 ( a2 ) 2 3 ( a2 3 ) 2 a1 ] + a 0 u 3 = p 3 e iϕ v 3 = p 3 e iϕ (3.8) D tan ϕ = D < 0 q D < 0 p u 1 = pe iϕ/3 v 1 = pe iϕ/3 u 2 = u 1 e 2πi/3 v 2 = v 1 e 2πi/3 (3.9) v 3 = v 1 e 2πi/3 u v (3.5)

3. 2 y 4 + a 3 y 3 + a 2 y 2 + a 1 y + a 0 = 0 y = x a 3 4 n p n (x) = a n x n + a n 1 x n 1 + + a 1 x + a 0 (3.13) p n (x) x k p n (µ) = ( ((a n µ + a n 1 )µ + a n 2 )µ a k x 4 + px 2 + + a 1 )µ + a 0 + qx + r = 0 (3.10) ( a3 ) 2 b n+1 = 0 p = a 2 6 4 b k = a k + µb k+1 (3.14) ( a3 ) ( a3 ) 3 q = a 1 2a 2 + 8 k = n, n 1,, 0 4 4 ( a3 ) ( a3 ) 2 ( a3 ) 4 b r = a 0 a 1 + a 2 3 0 4 4 4 p n (µ) (3.10) (x 2 + α) 2 (βx q 2β )2 = 0 p n (µ) = b 0 b k (3.10) p n (x) x µ p n (x) = (x µ)(b n x n 1 + b n 1 x n 2 + 2α β 2 = p α 2 q2 4β 2 = r (3.11) + b 2 x + b 1 ) + b 0 (3.15) β (3.13) (3.15) b k (3.14) 4(α 2 r)(2α p) = q 2 (3.14) b k p α n (x) x µ α β (3.11) α, β b n x n 1 + + b 2 x + b 1 x 2 + βx + α q 2β = 0 x 2 βx + α + q 2β = 0 (3.12) (3.10) p n (x) = d n (x µ) n + d n 1 (x µ) n 1 + + + d 0 x µ p n (x) p n (x) = (x µ)[(x µ)(b n x n 2 + n 1 d k x = w + µ w

3. 3 x k x ε k p n (x) = (x 2 ux v)(b n x n 2 + b n 1 x n 3 + + b 3 x + b 2 ) + b 1 (x u) + b 0 (3.16) y 1 = f(x 1 ) = f(x + ε 1 ) = f (x)ε 1 + 1 2 f (x)ε 2 1 b n+1 = b n+2 = 0 b k = a k + ub k+1 + vb k+2 (3.17) k = n, n 1,, 0 u = 2µ v = µ 2 p n (x) = (x µ) 2 (b n x n 2 + b n 1 x n 3 + + b 3 x + b 2 ) + b 1 (x 2µ) + b 0 p n (µ) = b 0 µb 1 p n(µ) = b 1 (3.18) x k = x + ε k (3.21) f(x) = 0 (3.20) x 3 ε 3 = f (x) 2f (x) ε 1ε 2 (3.22) x 2 ε 2 ε 1 f (x) 2f (x) ε 2 y = f(x) = 0 (3.19) x 1 x 3 x 2 x 1 x 3 x 2 f(x) (x 1, x 2 ) f(x) f(x) y 1 + y 2 y 1 x 2 x 1 (x x 1 ) f(x) = 0 x 3 = x 1 x 2 x 1 y 1 (3.20) y 2 y 1 (x 1 y 2 x 2 y 1 )/(y 2 y 1 ) y 3 = f(x 3 ) y 1 y 2 f(x) (robust ) x 4 = x 2 + x 3 2 y 3

3. 4 (3.22) ( ) f 2 ( f f(x) x f(x) ε 4 = 2f f f ) 3f (ε 1 ε 2 ) 2 (3.28) x = f 1 (y) y (x 1, x 2, x 3 ) 4 x 4 x [1234] (y 1, y 2, y 3 ) (3.26) x (y y 2)(y y 3 ) (y 1 y 2 )(y 1 y 3 ) x 1 + (y y 3)(y y 1 ) (y 2 y 3 )(y 2 y 1 ) x 2 3 (3.26) + (y y 1)(y y 2 ) (y 3 y 1 )(y 3 y 2 ) x 3 (3.23) y = y i x = x i y = f(x) = 0 x y = 0 k x k y 2 y 3 x 4 = (y 1 y 2 )(y 1 y 3 ) x y 3 y 1 1 + (y 2 y 3 )(y 2 y 1 ) x x 2 y 1 y 2 + (y 3 y 1 )(y 3 y 2 ) x 3 (3.24) f(x k + x) = 0 f(x k ) + f (x k ) x + = 0 x 1, x j x [1j] x 1 x j x 1 x y 1 j = 2, 3 (3.25) k+1 = x k + x = x k f(x k) f y j y (3.29) (x k ) 1 (3.20) x 3 x [12] x 3 x 4 x 4 = x [123] x [12] x [13] x [12] y 3 y 2 y 2 (3.26) k ε k ε k+1 = f (x) 2f (x) ε2 k (3.30) (3.24) (3.26) (3.26) (3.25) x 1, x 2, x 3 (3.24) x 4 0 ε 4 = f ( f 2f f f ) ε k+1 ε k 3f ε 1 ε 2 ε 3 (3.27) f(x) x x 1, x 2, x 3 x 3 x 1 x 2 x 3 = x f (x)f(x) < 1 [12] [f (x)] 2

3. 5 x f 2 x f 2 x = a 3 a 0 0 f 2 f(x) = x 2 a = 0 (3.29) x k+1 = x k x2 k a = 1 ) (x k + axk 2x k 2 Fortran C x = ϕ(x) (3.31) 2 3 10 x k k 2 3 10 0 1.0 1.5 5.0 1 1.5 1.75 3.5 x k+1 = ϕ(x k ) (3.32) 2 1.41666666 1.73214286 3.17857143 3 1.41421569 1.73205081 3.16231942 4 1.41421356 3.16227766 ϕ(x) ϕ(y) q x y 0 < q < 1 (3.33) a f 1 (x) = x 3 a = 0 f 2 (x) = x 2 a x u = nt + e sin u f 1 : x k+1 = 1 ( 2x k + a ) 3 x 2 k f 2 : x k+1 = x k(1 + 2a/x 3 k ) n t e 2 + a/x 3 k u 3 u t 10 ϕ(u) = nt + e sin u k f 1 f 2 e 1 0 3.0 3.0 ϕ(u) (3.33) 1 2.37037037 2.20312500 2 2.17350863 2.15445072 3 2.15460159 2.15443469 4 2.15443470 f 2 (x) = 2(x3 a) x 3 f 1 x k+1 x k q x k x k 1 (e = 0.0934) nt = π/6

3. 6 k u k 0 0.5235987 1 0.5702987 b 0 2 0.5488112 u u + b 0 v v = b 0(u, v) b 3 0.5479608 1 u u + b 1 v v = b 1(u, v) 4 0.5479269 5 0.5479256 b 0 / u c k = b k d k = b k 1 (3.35) u v (3.17) u c n+1 = c n+2 = 0 c k = b k + uc k+1 + vc k+2 (3.36) k = n, n 1,, 0 v d n+1 = d n+2 = 0 d k = b k + ud k+1 + vd k+2 n p n (x) x 2 ux v k = n, n 1,, 0 p n (x) = (x 2 ux v)(b n x n 2 + b n 1 x n 3 + + b 2 ) + b 1 (x u) + b 0 (3.34) b k (3.17) c 1 u + c 2 v = b 1 (3.37) b 1, b 0 0 x 2 ux v p n (x) u v p n (x) u, v 6 b 0 (u, v) = 0 b 1 (u, v) = 0 (u, v) +4x 3 + 3x 2 + 2x + 1 (3.38) (u, v) ( u, v) b 0 b 1 b 0 (u + u, v + v) = 0 b 1 (u + u, v + v) = 0 d k = c k c 0 u + c 1 v = b 0 p 6 (x) = 7x 6 + 6x 5 + 5x 4

3. 7 k u v b 0 b 1 0 1.00000e+0 2.00000e+0 5.14000e+2 2.55000e+2 1 8.17127e 1 5.12230e+0 1.46922e+3 3.31247e+2 2 6.30386e 1 3.54868e+0 4.55822e+2 9.81376e+1 3 5.28856e 1 2.36285e+0 1.38267e+2 2.84699e+1 4 4.78485e 1 1.51738e+0 4.15141e+1 7.84127e+0 10 1.22970e+0 4.14370e 1 6.34216e 1 5.71124e 1 11 1.16159e+0 4.79867e 1 1.88651e 3 1.71509e 2 12 1.26827e+0 4.84868e 1 2.24267e 4 4.07016e 4 13 1.26822e+0 4.84843e 1 1.83115e 8 3.86946e 8 n 2 p n (z i ) z i = a n j i (z (3.40) i z j ) i = 1, 2,, n (u, v) (3.29) (u, v) z i = p n(z i ) p n(z i ) (3.39) p n(z) z = z i (3.40) (3.40) z DKA n i p n (z) z i (i = 1, 2,, n) z i (3.39) p n (z) = a n (z z 1 z 1 )(z z 2 z 2 ) (z z n z n ) (3.39) p n (z) = z z j j z j p n (z) = a n a n n i=1 n i=1 (z z i ) z i (z z j ) j i p n (z i ) p n(z i ) z = z i (3.18) p n (z i ) 0 (3.40) p n (z i ) = z i a n (z i z j ) j i p n(z) z = z i 1 z i = j i 1 1 p n(z i ) z i z j p n (z i ) (3.41) z i

3. 8 P (re iϕ ) = r n e inϕ n 1 0 P (w) = 1 ( p n w a ) +c n 2 r n 2 e (n 2)iϕ + + c 0 = 0 n 1 a n na n e inϕ = w n + c n 2 w n 2 + + c 0 (3.42) r n + c n 2 r n 2 e 2iϕ + + c 0 e niϕ = 0 c k n 1 r > r 0 S(w) = w n c n 2 w n 2 c n 3 w n 3 c 1 w c 0 (3.43) 0 +( c n 3 + c n 3 e 3iϕ )r n 3 S(w) + + ( c 0 + c 0 e niϕ ) = δ w > 0 S(0) < 0 S( ) > 0 P (w) w > r 0 S(w) w = r 1, r 2 r n 2 S(r 2 ) r n 1 S(r 1 ) = c n 2 (r 2 1 r 2 2 ) + c 0 (r n 1 r n 2 ) = 0 r 1 r 2 S(r) > 0 r > r 0 0 r r S(r) w > 0 1 r 0 r 0 S(r) > 0 r S(r 0 ) = 0 r 0 > 0 (3.44) p n (z) r z j = a n 1 na n S(r) > 0 r > r 0 { [ 2(j 1)π +r exp i + π ]} (3.46) P (w) n 2n j = 1, 2,, n r 0 w r 0 (3.45) S(r) = r n c n 2 r n 2 c 0 δ > 0 ( c n 2 + c n 2 e 2iϕ )r n 2 (3.40) r 0 w = re iϕ r > r 0 (3.40) zi

3. 9 δ k (3.40) Durand Kerner (3.46) n Aberth x j M DKA 6 (3.38) (3.46) (3.41) p n (x) = a n (x x j ) j=1 4 x k n M M p n (x) M -1.0 0 1.0 n M M 3 3 x 1 < x 3 < x 2 f(x) x 1 < x < x 2 n p n (x) a k δ k f(x 3 ) < f(x 1 ) f(x 3 ) < f(x 2 ) x 0 x 4 x 3 (a n + δ n )x n + (a n 1 + δ n 1 )x n 1 + x 4 + (a 0 + δ 0 ) = 0 ( ) x = x 0 + ε ε, δ k f(x) p n (x 0 + ε) + n δ k x k 0 = 0 k=1 p n (x 0 ) + p n(x 0 )ε + k k ε = δ kx k 0 p n(x 0 ) δ k x k 0 = 0 (3.47) x 1 x 3 x 2

3. 10 x 3 x 3 x 4 x 3 x 1 x 2 x 1 = u x 4 x 3 x 2 x 1 = v x 1 x 4 x 3 x 2 (x 1, x 3, x 2 ) x 3 0.382 : 0.618 x 4 x 3 u + v = 1 u f(x 4 ) f(x 4 ) > f(x 3 ) x 3 x 4 (x 1, x 3, x 4 ) u = v 1 u 0.62 u 2 3u + 1 = 0 u = 3 5 = 0.381966 (3.48) 2 v = 5 2 = 0.236068 u v x 1 x 3 x 4 x 2 f(x 4 ) < f(x 3 ) (x 3, x 4, x 2 ) x 4 2 1 1.62 x 3 A B 0.382 : 0.618 = 1 : 1.618 1 2 1 1.29 x 3 1 1.44 1 1.65 x 3 v x 4

4. 1 4 a (3) 33 x 3 = b (3) 3 (f) a 11 x 1 + a 12 x 2 + a 13 x 3 = b 1 (a) a (2) 22 x 2 + a (2) 23 x 3 = b (2) 2 (d) a 11 x 1 + a 12 x 2 + a 13 x 3 = b 1 a 21 x 1 + a 22 x 2 + a 23 x 3 = b 2 a 31 x 1 + a 32 x 2 + a 33 x 3 = b 3 (a) (b) (c) (a) (b) x 1 (a) a 21 /a 11 (b) a (2) 22 x 2 + a (2) 23 x 3 = b (2) 2 (d) l 21 = a 21 a 11 a (2) 22 = a 22 l 21 a 12 a (2) 23 = a 23 l 21 a 13 b (2) 2 = b 2 l 21 b 1 (c) x 1 x 1, x 2 x k 1 (b) (c) (a) (a) a 31 /a 11 a (1) 11 x 1 + a (1) 12 x 2 + a (1) 13 x 3 + + a (1) 1n x n = b (1) 1 (c) l 31 = a 31 a 11 a (2) 32 = a 32 l 31 a 12 a (2) 33 = a 33 l 31 a 13 b (2) 3 = b 3 l 31 b 1 a (2) 32 x 2 + a (2) 33 x 3 = b (2) 3 (e) x 1 x k (d) (e) x 2 l 32 = a(2) 23 a (2) 22 a (3) 33 = a(2) 33 l(2) 32 a(2) 23 b (3) 3 = b (2) 3 l (2) 32 b(2) 2 a (3) 33 x 3 = b (3) 3 (f) (a) (d) (f) x 3 = b(3) 3 a (3) 33 x 2 = 1 a (2) 22 ( ) b (2) 2 a (2) 23 x 3 x 1 = 1 a 11 ( b 1 a 12 x 2 a 13 x 3 ) n (g) a i1 x 1 + a i2 x 2 + + a in x n = b i (4.1) i = 1, 2,, n a (2) 22 x 2 + a (2) 23 x 3 + + a (2) 2n x n = b (2) 2.... (4.2) a (k) kk x k + + a (k) kn x n = b (k) k a (k) k+1k x k + + a (k) k+1n x n = b (k) k+1.... a (k) nk x k + + a (k) nn x n = b (k) n a (1) ij = a ij b (1) i = b i k k + 1 n l ik = a(k) ik a (k) kk a (k+1) ij b (k+1) i = a (k) ij = b (k) i k < i j n (4.3) l ik a (k) kj (4.4) l ik b (k) k (4.5) 0 k < n

4. 2 Fortran do k=1, n-1 do i=k+1, n lik=a(i,k)/a(k,k) do j=k+1, n a(i,j)=a(i,j)-lik*a(k,j) enddo b(i)=b(i)-lik*b(k) enddo enddo (4.4) k (n k) 2 a (1) 11 x 1 + a (1) 12 x 2 + a (1) 13 x 3 + + a (1) 1n x n = b (1) 1 ( ) a (2) 22 x 2 + a (2) 23 x 2 + + a (2) 2n x n = b (2) 2 x n = b(n) n a (n) nn x i = 1 a (i) ii a (3) 33 x 3 + + a (3) 3n x n = b (3) 3 (4.6)..... ( b (i) i n j=i+1 a (n) nn x n = b (n) n ) a (i) ij x j (4.7) i = n 1, n 2,, 1 (4.1) (4.4) a (k) ij 1 3 n3 + O(n 2 ) (4.8) a (k+1) ij a (k) ij a (k) n ij a ij n = 100 n = 1000 (4.5) b (k) i b i (4.3) a (k) ik (4.8) l ik 1/3 2/3 a ij b i (4.5) (4.7) 1 0 n (4.4) (4.5) (4.7) n Fortran (4.4) k = 1 n 1 n (4.4) n 1 (n k) 2 1 3 n3 (4.3) kk ( ) k 0 k=1 (4.3) (4.5) (4.7) n (4.3), (4.5), (4.7) 1 2 n2 a (k)

4. 3 0 x 3 x 3 x 1 x 2 x 1 x 2 x 1 x 2 (4.10) k k a (k) ik k i n (4.10) i i k x i 0 LU a (k) k ij n n A (k) A (k) A (k+1) a (k) ij k i, j n A (k+1) = M k A (k) i, j i k j k M k n n k = 1 M 1 1 0 0 l 21 1 0 M 1 = det[a ij ] = a (1) 11 a(2) 22 a(k) kk... a(n) nn (4.9)...... l n1 0 1 M k (k, k) [ l k+1,k, l k+2,k,, l n,k ] T ( T ) n 1 0 A U n = 3 x 1 a (2) 22 = a(2) 32 = 0 M n 1 M 2 M 1 A = U (4.11) U (4.6) a (2) 23 x 3 = b (2) 2 a (2) 33 x 3 = b (2) 3 b (2) 2 a (2) 23 = b(2) 3 a (2) 33 (4.10) L = M 1 1 M 1 2 M 1 n 1 (4.12) (4.11) A A = LU (4.13)

4. 4 A LU L x 2 x 3 M 1 k M k l i,k 1 0 0 l 21 1 0 L = l 31 l 32 1 0........... l n1 l n2 l n3 1 (4.14) L A A (4.13) LU l ij n (4.3) Fortran A det A = det L det U det L = 1 A (4.9) a(i,j) det=1 do k=1, n akk=a(k,k) det=det*akk a(k,k)=1 n = 3 do j=1, n+1 a 11 a(k,j)=a(k,j)/akk x 1 + a 12 x 2 + a 13 x 3 = b enddo 1 do i=1, n if( i.ne.k ) then aik=a(i,k) a 21 a 31 a(i,k)=0 x 1 do j=1, n+1 a 22 a(i,j)=a(i,j)-aik*a(k,j) enddo x 1 + a 12 x 2 + a 13 x 3 = b 1 x 2 + a 23 x 3 = b 2 a 32 x 2 + a 33 x 3 = b 3 a 12 x 1 = b 1 x 2 = b 2 x 3 = b 3 n + 1 b i a(i,n+1) a 32 endif enddo enddo x 2

4. 5 LU (4.6) 0 A A = L U u 11 u 12 u 13 u 1n 0 u 22 u 23 u 2n U = 0 0 u 33 u 3n....... 0 u nn (4.15) L (4.14) L U L U = A (1,1) γ i x i 1 + α i x i + β i x i+1 = b i (4.18) i = 1, 2, 3,, n γ 1 = β n = 0 α i, β i, γ i a ij a i x i + β i x i+1 = c i i = 1, 2,, n (4.19) u 11 = a 11 a i, c i a 1 = α 1 c 1 = β 1 l i1 = a a i = α i γ i β i 1 (4.20) i1 a i 1 i = 2, 3,, n u 11 c i = b i γ i c i 1 a L U i 1 i = 2, 3,, n (1,2) (2,2) u 12 = a 12 u 22 = a 22 l 21 u 12 l 21 x i = 1 (c i β i x i+1 ) (4.21) u 12 u 22 a i i = n, n 1,, 1 l i2 = 1 ( ) a i2 l i1 u 12 i = 3, 4,, n u 22 k k + 1 j = 1, 2,, n i 1 l ij = a ij l ik u kj i = 1, 2,, j (4.16) k=1 l ij = 1 ( a ij u jj j 1 ) l ik u kj k=1 i = j+1, j+2,, n (4.17) a i, c i (4.12) x i α i > β i + γ i (4.22) LU (4.1) n n A = [a ij ] b i b = [b 1, b 2,, b n ] T ( T ) x = [x 1, x 2,, x n ] T (4.1) ( ) Ax = b (4.23)

4. 6 (4.31) (4.29) A 1 x = A 1 b (4.24) x = 0 (4.32) x x x 0 x = 0 x = 0 (4.25) αx = α x x + y x + y A λ i α ρ(a) = max λ i (4.34) L 2 ( n x 2 = x i 2) 1/2 i=1 x 1 = L 2 (4.26) n x i L 1 (4.27) i=1 x = max 1 i n x i L (4.28) A b Ax A = sup x 0 x (4.29) x Ax = A x b (4.29) L 1 L A 1 = max a ij j A = max i A(x + δx) = b + δb δx = A 1 δb (4.25) i a ij (4.33) j L 2 L 1 L L 2 i A A L 2 A 2 = ρ(a) (4.35) Ax = b b δb x + δx AB A B (4.30) Ax A x (4.31) b = Ax A x δx = A 1 δb A 1 δb

4. 7 δx x κ(a) δb b (4.36) κ(a) = A A 1 (4.37) κ(a) κ(a) A κ (4.1) A 1 1 δa x + δx (A + δa)(x + δx) = b δx = A 1 δa(x + δx) U T z = e L T v = z (4.40) δx x + δx κ(a) δa A (4.38) e v δa / A A L 2 e 1 A 1 A 1 (4.35) κ 2 (A) = max i λ i min i λ i (4.39) (3.29) A 1 1 = A T = sup A T e e A 1 A T A T e = v e = max e i = 1 (4.40) A A z 1 = e 1 = 1 u 11 u 11 ( ) z k = 1 k 1 e k u ik z i (4.41) L 2 k = 2, 3,, n A A 1 e k = ±1 z k L 1 A 1 (4.33) sign(x) A 1 ( ) 1 k 1 A 1 e k = sign u ik v i z (4.40) (4.33) A A 1 = A T v n = z n i v = max i v i A LU v U T z z v e T = [1, ±1, ±1, ] 1 u kk i=1 i=1

4. 8 v k = z k n i=k+1 l ik z i (4.42) k = n 1, n 2,, 1 v A 1 1. = max v i (4.43) i κ 1 (A) = A 1 A 1 1 Ax = b x δx x + δx v A(x + δx) = b A 1 1 δx 1 2 3 4 5 6 2 2 3 4 5 6 3 3 3 4 5 6 A = 4 4 4 4 5 6 5 5 5 5 5 6 6 6 6 6 6 6 Aδx = b Ax = r (4.44) x r x + δx L 1 (4.43) n n 3 /3 κ 1 144 x LU 27.72 (4.43) 0.2682 L 2 n 2 κ 2 = 103 A n < 10 n A

5. 1 5 x 0, x 1,, x n f[x] = f(x) f(x) f n (x) f[x 0, x 1,, x k, x] = f(x) f[x 0, x 1,, x k 1, x] f[x 0, x 1,, x k 1, x k ] f n (x k ) = f(x k ) k = 0, 1,, n (5.1) x x k (5.3) f n (x) x a k = f[x 0, x 1,, x k ] 3 Π n (x) = (x x 0 )(x x 1 ) (x x n ) n = (x x i ) i=0 n Π n (x) f n (x) = (x x k )Π k=0 n(x k ) f(x k) (5.2) (5.1) (5.3) (5.4) 1 f(x 0 ) x 0 x (5.2) I 0i (x) = x i x 0 f(x i ) x i x (5.2) i = 1, 2,, n f n (x) x n 1 I 01 (x) x 1 x I 01i (x) = f n (x) = a 0 + a 1 (x x 0 ) x i x 1 I 0i (x) x i x +a 2 (x x 0 )(x x 1 ) + i = 2, 3,, n f n [x 0, x] f n(x) f n (x 0 ) x x 0 f n [x 0, x] = a 1 + a 2 (x x 1 ) + f n [x 0, x 1 ] = a 1 +a n (x x 1 ) (x x n ) f n (x) = f[x 0 ] + f[x 0, x 1 ](x x 0 ) +f[x 0, x 1, x 2 ](x x 0 )(x x 1 ) + +f[x 0, x 1,, x n ](x x 0 )(x x 1 ) (x x n ) +a n (x x 0 )(x x 1 ) (x x n ) 1 I 012 k 1i (x) = x i x k 1 I 012 k 1 (x) x k 1 x f n (x 0 ) = a 0 I 012 k 2i (x) x i x i = k, k+1,, n (5.4) (5.5) x(i) x i y(i) f(x i ) x = x c y c real do x(0:n), y(0:n), f(0:n) i=0, n f(i)=y(i) enddo do k=0, n-1

5. 2 do i=k+1, n f(i)=f(k)+(f(k)-f(i)) * *((x(k)-xc)/(x(i)-x(k))) enddo enddo yc=f(n) f(i) (5.2) + x j x f j 1 + x x j 1 f j 1 j n (5.6) (5.4) (5.5) h j h j x n x j 1 < x < x j h j = x j x j 1 f j = f(x j ) x 0, x 1,, x n S j (x j ) = f j S j (x j 1 ) = f j 1 [ 1 f(x) = S j(x) = h ( ) 2 j 1 + 25x 2 6 M xj x j 1 3 + 1] h j [ x = 1.0, 0.8,, 0.8, 1.0 11 + h ( ) 2 j x 6 M xj 1 j 3 1] + f j f j 1 S x j x x x j 1 j (x) = M j 1 + M j h j h j 10 M j x j 1.0 x j 1 x x j [ (xj ) 3 ( ) ] S j (x) = h2 j 6 M x xj x j 1 h j h j [ (x ) 3 ( ) ] + h2 j 6 M xj 1 x xj 1 j h j h j x j S j(x j ) = S j+1(x j ) h j h j 0.5-1.0-0.5 0 0.5 1.0 h j 6 M j 1 + h j + h j+1 M j + h j+1 3 6 M j+1 = f j+1 f j f j f j 1 (5.7) h j+1 h j M 0 = M n = 0 j = 1, 2,, n 1 M 1 M n 1 4 x 0 < x 1 < (5.7) < x n

5. 3 a i n + 1 0 a i (5.7) M j S n = f 2 (5.6) f(x) [a, b] ϕ i (x) u(x) v(x) L 2 (5.11) ϕ i (x) (u, v) = b a u(x)v(x)dx (5.8) u 2 a 2 = (u, u) (5.9) i n u 2 (4.25) n = 10 n = 11 L 2 (5.11) 2 { 1 i = j (ϕ i, ϕ j ) = δ ij = (5.12) 0 i j ϕ i (x) f(x) [a, b] δ ij ϕ i (x) n f(x) n f n (x) = a i ϕ i (x) (5.10) i=0 a i = (f, ϕ i ) (5.13) b n S n = f(x) a i ϕ i (x) 2 dx n a i=0 n S n = f 2 a 2 i a k i=0 b n [f(x) a i ϕ i (x)]ϕ k (x)dx = 0 n S n ( a i=1 ) n a i n a i (f, ϕ i ) = f 2 f n 2 lim f f n = 0 (5.14) n n a i (ϕ i, ϕ k ) = (f, ϕ k ) (5.11) i=0 f(x) (f, ϕ i )ϕ i (x) (5.15) k = 0, 1,, n i=0 x i (5.11) (5.10) (5.11) i = k i=0

5. 4 f(x) (f, ϕ i ) (5.15) x f(x) n = 5 1 ϕ i cos kx, sin kx [ π, π] π f(x) 0 π f(x) 1 2 a 0 + a k = 1 π b k = 1 π (a k cos kx + b k sin kx) (5.16) k=1 π π π π f(x) cos kxdx f(x) sin kxdx 1/π cos kx, sin kx f(x) x π f n (x) b k f n (x + π) = f n (x) f n (x) = 1 2 + n [ 1 π ] f(y) sin kydy sin kx π k=1 0 f n (x) x f(x) = 1 2 + 1 π n [cos k(x y) 2π 0 k=1 cos k(x + y)]dy 1 0 < x π x f(x) = 1/2 x = 0 (5.17) x y 0 π x < 0 df n (x) = 1 π n d [cos k(x y) dx 2π 0 dy k=1 + cos k(x + y)]dy 2 a 0 = 1 b 2k 1 = (2k 1)π df n (x) = 1 n [cos kx(1 cos kπ)] dx π a 0 a k b k 0 k=1 df 2n 1 (x)/dx f 2n 1 (x) = 1 2 + 2 π n k=1 sin(2k 1)x 2k 1 n = 10 sin 2nx n = 5, 10 = π sin x π df 2n 1 (x) dx = 2 π 1 0 n cos(2k 1)x k=1 π (5.18)

5. 5 df 2n 1 (x)/dx 0 x = mπ 2n m = ±1, ±2,, 2n 1 f 2n 1 (x) x = 0 n sin x x = π/2n Sinc(x) =. 1 x2 n 2 3! + x4 5! f 2n 1 (x) (5.18) x 6 /7! f 2n 1 (π/2n) = 1 2 + 1 π = 1 2 + 1 2nπ π 0 π/2n sin 2ny 0 sin y dy sin y sin(y/2n) dy n 1 π sin y 2nπ sin(y/2n) dy 1 π sin y π y dy 0 0 (Lanczos, Si(x) = x 0 sin y y dy x = π n 1 Sinc(x) =. 1 0.16670387x 2 + x4 (5.21) 120 (5.20) (5.21) x 2 Sinc(x) = sin x x (5.19) x = 0 0/0 sin x x = 0 sin x sin x/x x = 0 sin x (5.20) x = 0 (5.20) x = π/2n x x 2 < 1/2 1 1 = 2.5 10 5 7! 23 0 ) x 2 1/2 x 2 = 1 2 sin θ Si(π) = 1.8519370 f 2n 1 (π/2n) x 6 7! = 1 7!2 3 sin3 θ = 1 ( 3 7!2 3 4 sin θ 1 ) sin 3θ 4 f 2n 1 (π/2n) 1.089489 n sin θ sin 3θ 1 1 3 7!2 3 4 sin θ = 1 3 7!2 3 2 x2 x (5.20)

5. 6 3.0e-5 2.0e-5 1.0e-5-0.5 0 0.5 c k (5.24) x cos nθ cos nθ x n T n (x) = cos(n cos 1 x) (5.26) T 3 (x) = 4x 3 3x T 4 (x) = 8x 4 8x 2 + 1 sin nθ x n T n+1 (x) 2xT n (x) + T n 1 (x) = 0 (5.27) x 1 cos f(x) (5.24) x x = cos θ (5.22) cos f(cos θ) θ (5.16) (5.23) c k cos kθ f(cos θ) 1 x 1 f(x) 2 c 0 + c k cos kθ (5.23) k=1 f(x) = a n x n π c k = 1 f(cos θ) cos kθdθ π π a n x cos kθ 1 c k n x n = 1 [n/2] ( ) n 2 n 1 T n 2k (x) (5.28) k k=0 f n (x) = 1 n 2 c [n/2] n/2 0 + c k cos kθ (5.24) n T 0 (x) k=1 1/2 ( n + 1 p.89 ) c n+1 (5.24) f(x) = c k T k (x) x 8 c k c k V n = V n+1 = 0 V k 1 = c k + 2xV k V k 1 (5.25) k = n, n 1,, 1 f n (x) = 1 2 c 0 + xv 0 V 1 x n T 0 (x) = 1 T 1 (x) = x T 2 (x) = 2x 2 1 n=0 k=0 (5.27)

5. 7 Sinc Sinc (5.19) θ k n n cos mθ k = sin mθ k = 0 m 0 (5.30) (5.19) x π θ (5.20) k x T n (x) 1 x 1 n Sinc(π x) = a n x n n=0 a n = ( 1)n π 2n (2n + 1)! a n c n n 0 i j T i (x k )T j (x k ) = n/2 i = j 0 (5.32) n a n c n k=1 n i = j = 0 0 1.000000e+0 1.415723e+0 i, j < n 1 1.644934e+0 1.789468e+0 T i (x) T j (x) x k 2 8.117425e 1 4.190149e 1 3 1.907519e 1 4.842341e 2 1 x 1 f(x) 4 2.614786e 2 3.296368e 3 5 2.346082e 3 1.473941d 4 c k = 2 n f(x j )T k (x j ) (5.33) 6 1.484289e 4 4.654237e 6 n 7 6.975877e 6 1.092549e 7 8 2.531219e 7 1.980867e 9 9 7.304716e 9 2.85340e 11 f n (x) f n (x k ) = f(x k ) k=1 Sinc (5.33) (5.34) (5.32) ( ) [a, b] f(x) n Sinc p n (x) = a 0 + a 1 x + + a n x n c k cos nθ 0 θ π ε = max p n(x) f(x) (5.35) a x b n θ k = (k 1/2) π (a 0, a 1,, a n ) k = 1, 2,, n (5.29) n k=1 T n (x k ) = 0 x k = cos θ k (5.31) (5.30) j=1 f n (x) = 1 n 1 2 c 0 + c k T k (x) (5.34) k=1

5. 8 tan x/x n + 2 ( ) x = π y 1 y 1 4 y 4 tan x/x (5.21) ( y tan tanh T 5 (y) ) 4 (5.36) y = 1 n T n+1 (x) T n+1 (x) x 1.0e-4 x 1 f(x) n T n+1 (x) -1.0 0 1.0 y T n+1 (x) (5.29) n n + 1 n + 1 0 0-1.0e-4 a i x j a 0 + a 1 x j + a 2 x 2 j + + a n x n = f(x j ) (5.36) tan x/x j = 1, 2,, n + 1 tan x = 1 + 1 x 3 x2 + 2 15 x4 + 17 315 x6 n + 1 a i + 62 2835 x8 + (5.38) x 8 e(x) = p n (x) f(x) y 0.3 e(x) x k y 4 e(x k ) x k n + 2 e(x k ) e(x k ) (5.8) ε w(x) 0 p n (x k ) = f(x k ) ± ε (5.37) k = 1, 2,, n + 2 p k (x) a i ± e(x k ) (p j, p k ) = w(x)p j (x)p k (x)dx = λ k δ jk a (5.39) p k (x) x k µ k 0 µ 0 = 1 b

5. 9 [a, b] w(x) λ n α k β k γ k P n (x) [ 1, 1] 1 2/(2n + 1) 2 1/k 0 1 1/k T n (x) [ 1, 1] 1/ 1 x 2 π/2 2 0 1 L n (x) [0, ) exp( x) 1 1/k 2 1/k 1 1/k H n (x) (, ) exp( x 2 ) π2 n n! 1 0 k 1 p k (x) (5.41) x = x j y x j y = x l p 0 (x) = µ 0 p 1 (x) = 0 (5.32) p k (x) = µ k x k + p k (x) = (α k x β k )p k 1 (x) γ k p k 2 (x) α k = µ k β k = α k (xp k 1, p k 1 ) (5.40) µ k 1 λ k 1 γ k = α c k = 1 n w j p k (x j )f(x j ) kλ k 1 λ k j=1 α k 1 λ k 2 (5.43) α k β k n 1 p k (x) p k 1 (x) f n (x) = c k p k (x) (5.44) k=0 p k 2 (x) n 1 p k (x)p k (y) = µ n 1 λ k µ n λ n 1 (5.43) k=0 p n(x)p n 1 (y) p n 1 (x)p n (y) [ ] (5.41) n n 1 p k (x j )p k (x) x y f n (x) = w j f(x j ) x = x l (5.42) x j = l f n (x l ) = f(x l ) x n (5.44) p n (x), p n 1 (x), p k (x) (5.40) (5.41) n p n (x) f n (x) = p k (x) p n (x) (x x j=1 j )p n(x j ) f(x j) (5.45) x 1, x 2,, x n (5.41) n 1 1 [p k (x j )] 2 w j λ k k=0 n 1 k=0 = µ n 1 µ n λ n 1 p n 1 (x j )p n(x j ) (5.42) p k (x j )p k (x l ) λ k = 0 j l j=1 k=0 λ k (5.2) x j p n (x) ( )

5. 10 f(x) f(x) = b 0 + b 1 + a 1 a 2 b 2 + a 3 b 3 +... (5.46) N+1 = b N+1 a i, b i x n = N, N 1,, 0 0 f(x) f(x) = b 0 + a 1 a 2 a 3 a 4 b 1 + b 2 + b 3 + b 4 + (5.47) N (5.47) tan x = x x 2 x 2 1 3 5 e x = 1 + x x x x 1 2+ 3 2+ x 2n + 1 log(1 + x) = x x x 2x 2x 1+ 2+ 3+ 2+ 5+ nx nx 2+ 2n + 1+ tan 1 x = x x 2 4x 2 1+ 3+ 5+ n 2 x 2 2n + 1+ x 2+ f n (5.47) a n, b n j = 1, 2 Z ν 1 (x) + Z ν+1 (x) = 2ν x Z ν(x) f j f j 1 Z ν 1 (x) Z ν (x) Z ν+1 (x) Z ν (x) = 2ν x 1 Z ν (x)/z ν+1 (x) = 2ν x 1 1 2(ν + 1)/x 2(ν + 2)/x 1 2(ν + 3)/x = x 2(ν + 1) x 2 2(ν + 2) n+1 n = b n + a n+1 n+1 (5.48) N A 1 = 1 B 1 = 0 A 0 = b 0 B 0 = 1 A j = b j A j 1 + a j A j 2 (5.49) B j = b j B j 1 + a j B j 2 j = 1, 2, f n = A n B n (5.49) A j, B j C j = x 2 2(ν + 3) f j = f j 1 C j D j C j D j A j A j 1 D j = B j 1 B j (5.47) f 0 = C 0 = b 0 D 0 = 0 n C j = b j + a j (5.50) C j 1 n = b n + a n+1 1 D j = j = 1, 2, b n+1 + b j + a j D j 1

5. 11 C j D j 1 f(x) m n C j D j 0 m f(x) A a i x i m(x) B n (x) = i=0 (5.51) n 1 + b j x j j=1 (5.49) a i b j (5.50) a i b j a i b j b 0 = 1 (5.50) m+n+1 (5.49) x i a 0 + a 1 x i + a 2 x 2 i + + a m x m i tan x 3 tan x x. = 1 1 x 2 3 x 2 5 = 15 x2 15 6x 2 a 0 a m, b 1 b n x i tan x/x x tan x/x m + n + 1 0.10 1.003346 0.640e-9 0.009e-9 0.20 1.013550 4.204e-8 0.923e-9 0.30 1.031121 5.001e-7 5.432e-8 0.40 1.056983 2.990e-6 9.938e-7 0.50 1.092605 1.239e-5 9.631e-6 0.60 1.140228 4.110e-5 6.275e-5 0.70 1.203269 1.182e-4 3.124e-4 e(x) = A m(x) B n (x) f(x) 0.80 1.287048 3.098e-4 1.285e-3 0.90 1.400176 7.675e-4 4.601e-3 m + n + 2 1.00 1.557408 1.852e-3 1.490e-2 x j tan x (5.38) x > 0.5 A m (x j ) = [f(x j ) ± ε] B n (x j ) (5.53) j = 1, 2,, m + n + 2 1 ε ± e(x) m x j ε n a i b j m + n + 2 e(x) = f(x i )(1 + b 1 x i + b 2 x 2 i + + b n x n i ) i = 1, 2,, m + n + 1 (5.52) x i

5. 12 (5.53) m+n+1 a i b j m + n + 2 (5.50) ε a i (5.52) tan x/x m = n = 2 A 2 (x) = 1 0.11111233x 2 + 1.0582954e-3x 4 B 2 (x) = 1 0.44444567x 2 + 1.5873516e-2x 4 A m (x), B n (x) f(x) 0 x 1 10 7 x 8 f(x) = c k x k x = 1 1.5 10 2 k=0 (5.51) tan x B n (x) c k x k A m (x) k=0 x = 2 x i a i, b j 10 4 n i b j c k x j+k = b j c i j x i tan x j=0 k=0 i=0 j=0 x m A 5 1 1 9 x2 + 1 945 x4 B i 5 1 4 9 x2 + 1 63 x4 b j c i j = a i i = 0, 1,, m (5.54) j=0 A m (x) B n (x) x m A m (x) B n (x) 0 A m (x)/b n (x) i c k b j c i j = 0 (5.55) j=0 c k (5.55) b 1 b n π/2 8 x > 1 tan x x = π/2 (5.49) j = 5 i = m + 1, m + 2,, m + n i = m + n b j b 0 = 1 b 1 b n n n i b j e b n min(i, n) x 1 + x/2 + x2 /10 + x 3 /120 1 x/2 + x 2 /10 x 3 /120

5. 13 1 + 10x + 5x 2 cos x 5 + 10x + x 2 3 5 + 35x + 14x 2 2T j (x)t k (x) = T j+k (x) + T j k (x) x 14 + 35x + 5x 2 tan 1 x 945 + 735x2 + 202x 4 x 945 + 1050x 2 + 225x 4 x 3 x x = 0 T j (x) c k T k (x) x ε 0.8 < x < 1.2 ε = 8 10 7 0.5 < x < 2.0 ε = 3 10 4 1 x 1 n b j c (j) k = 0 k = m+1, m+2,, m+n f(x) = c k T k (x) j=0 k=0 n b j c (j) k = a k k = 0, 1, 2,, m (5.56) (5.51) A m (x), B n (x) b j m n A m (x) = a i T i (x) B n (x) = b j T j (x) a i i=0 j=0 B n (x)f(x) 1963 : k=0 = 1 2 = k=0 c (j) k c k [T j+k (x) + T j k (x)] k=0 c (j) k T k(x) c (0) k = c k n B n (x)f(x) = b j c (j) k T k(x) j=0 j=0 k=0 m a i T i (x) i=0

6. 1 6 A = QR Q = [q 1, q 2,, q m ] (6.3) 4 R = [r ij ] A LU QR Q n m R m m r 11 r 12 r 1m r 22 r 2m QR R =.... 0 r mm q k Q n m (n m) A n a j A = [a 1, a 2,, a m ] (6.1) A n m QQ T a j (6.4) q 1 = a 1 r 11 = a 1 r 11 L 2 a 2 q 1 q 1 a 2 q j a 2 q 1 q k u 2 = a 2 q 1 r 12 r 12 = q T 1 a 2 a j T u 2 q 1 u 2 1 q 2 = u 2 r 22 j=1 r 22 = u 2 q k = u k r kk r kk = u k (6.2) u k a k r k 1,j = q T k 1a (k 1) j j = k, k+1,, m (6.5) q 1, q 2,, q k 1 q k = a(k) k r kk = a (k) k j=1 j=1 Q T Q = I m (6.4) I m m m A = [a (1) 1, a(1) 2,, a(1) m ] q 1 q 1 q 2 a (1) 2, a(1) 3,, a(1) m k q k q 1 a (2) 2 k 1 u k = a k q j r jk r jk = q T q k 1 j a k q k a (k) j = a (k 1) j r kk q k 1 r k 1,j u k q m k 1 k a k = q k r kk + q j r jk = q j r Q R (6.3) jk (6.4)

6. 2 (6.5) r kk a w a (k) P a = a w(w j (j k) T a) = [s, 0,, 0] T (6.9) s 2 = a 2 a (k) j a (k) k P a 0 a T Q A R a 2 (w T a) 2 = sa 1 Q a n = m Ax = b 1 a w T a Q T (6.9) w w Q T QRx = Q T b Rx = Q T b w T a = s(s a 1 ) s = a sign ( a 1 ) r jk w = 1 w T a [a 1 s, a 2,, a n ] T (6.10) A b a (1) n+1 (6.5) s s a 1 n (n + 1) R r i,n+1 n P = I + uut (6.11) n s(a 1 s) a (n+1) n+1 0 u T = [a 1 s, a 2,, a n ] b P P b = b + βu β = ut b s(a 1 s) QR 2 n m (n m) A w a A a 1 P P = I ww T w 2 = 2 (6.6) P P T = P P T P = I (6.7) P 1 P 1 A 1 2 0 P 1 A n 1 a P 2 P 1 A I P 2 a (1, 1) 1 0 P k (k, k) P k P a 2 = a T P T P a = a 2 (6.8) I k 1 0 P (6.6) P k = 0 I n k+1 ww T QR

6. 3 I k k k A = QL P k Q Q 0 A m(m n) [ ] R P m P m 1 P 1 A = 0 (6.12) R m m 0 a m 1 q m (n m) m m = n u m 1 = a m 1 l m,m 1 q m P m l m,m 1 = q T ma m 1 q m 1 = u m 1 l m 1,m 1 = u m 1 Q = P 1 P 2 P m (6.13) l m 1,l 1 Q (6.3) Q n n Q (6.12) m [ ] u k = a k q j l jk l jk = q T j a k (6.14) Q T R j=k+1 A = 0 q k = u k l kk = u k QQ T = I [ ] (6.14) A = QL R A = Q 0 n > m Q m QL (6.3) (6.12) m 1 R Q n 1 0 R A 0 (6.12) [ ] 0 Ax = b (n = m) A P k P m P m 1 P 1 A = L b P k P n 1 P n 2 P 1 Ax = Rx = P n 1 P n 2 P 1 b Q n n QL A Q R A = QR Q l 11 0 l 21 l 22 L =....... l m1 l m2 l mm q m = a m l mm r kk l mm = a m A Ax = λx (6.15) λ A x λ det(λi A) = 0

6. 4 P a T 1 B = P 1 AP (6.16) A B a 1 2 (w T a 1 ) 2 = a 2 (6.15) P 1 11 + sa 21 a 1 2 = a 2 11 + s 2 P 1 AP P 1 x = λp 1 x By = λy x = P y y = P 1 x w 0 w = [0,,,, ] T a 11 a 21 a 1 s = sign( a 21 ) n B w T a 1 = s(s a 21 ) (6.17) P A B w A B y w = 1 w T [0, a 21 s, a 31,, a n1 ] T (6.18) a 1 B = P 1 AP 1 B P (1, 1) 1 b 11 = a 11 b 12 = b 21 = s 0 P 0 A 0 B = P 1 AP = P AP A n n n 2 0 0 A B = 0. B = P T AP...... 0 α 1 β 1 0 0 β 1 α 2 β 2 0 A a 1. 0 β 2 α 2 β.. 3 = P 1. 0........ 0 P 1 a 1...... αn 1 β n 1 P 1 a 1 = a 1 w(w T a 1 ) = [a 11, s, 0,, 0] T 0 0 β n 2 α n i=2 a 2 i1 w s a 21 P 1 (6.11) P 1 = I + uut s(a 21 s) P = I ww T w 2 = 2 u T = [0, a 21 s, a 31,, a n1 ] QR P 1

6. 5 A P = P 1 P 2 P n 2 (6.19) 0 0 0 0 P T AQ =... 0.... (6.21).......... 0 0 P 1 A P = P 1 P 2 P n 1 1 0 Q = Q P 1 1 1 Q 2 Q n 1 0 0 (6.19) B = 0 (6.20).......... 0 0 A P 1 0 0 0 1 QR..... c s P =....... A 1 0 s c. P 1 (6.10).. 0 0 0 1 B = P 1 A (6.22) B b T 1 Q 1 = [b 11, s, 0,, 0] b T 1 Q 1 (6.17) (6.18) P T P = I B 0 0 0 P 1 AQ 1 =........... 0 p kk = p ll = c = cos ϕ p kl = p lk = s = sin ϕ B = P 1 AP = P T AP B = P T A B = B P

6. 6 B k, l B k, l ϕ B B cos ϕ sin ϕ { b kj = a kj cos ϕ a lj sin ϕ { b lj = a kj sin ϕ + a lj cos ϕ b ik = b ik cos ϕ b il sin ϕ b il = b ik sin ϕ + b il cos ϕ 1 j n 1 i n (6.23) b ik b il i = k, l A 1 B (k, l) c = s = ct 1 + t 2 b kl = b lk = (a kk a ll ) sin ϕ cos ϕ +a kl (cos 2 ϕ sin 2 ϕ) (6.24) a kl 0 0 ϕ (6.25) (k, l) (l, k) B 0 τ = a ll a kk A D ± t 0 s 0 0 τ > 0 A B = P T 1 AP t = τ + τ > 0 τ 2 + 1 b 2 ij = τ < 0 a 2 ij sign( τ) i,j i,j t = τ (1 + 1 + 1/τ 2 ) (6.23) sign b 2 kk + 2b 2 kl + b 2 ll = a 2 kk + 2a 2 kl + a 2 ll b kl 0 ϕ b 2 kk + b 2 ll > a 2 kk + a 2 ll a kl 0 b ll a kk ϕ (6.24) b kl = 0 2a kl tan 2ϕ = a ll a kk b kl = 0 t 2 + a ll a kk a kl t 1 = 0 t = tan ϕ t = a kk a ll 2a kl ± (akk ) 2 a ll + 1 (6.25) 2a kl P (6.23) (6.25) 2a kl (6.23) b kk = c 2 a kk + s 2 a ll 2csa kl c s a ll a kk = c2 s 2 cs a kl a ll b kk = a kk ta kl b ll = a ll + ta kl (6.26)

6. 7 (6.23) D λ j b kj = a kj s(a lj + γa kj ) P j p j b lj = a lj + s(a kj γa lj ) (6.27) γ = s j k, l Ap j = p j λ j 1 + c 0 a kl (k, l) n 2 /2 (k, l) = (6.19) (1, 2), (1, 3) (1, n), (2, 3), (2, 4) p n (λ) det(λi B) 0 λ n A k k (6.23) P 1, P 2, p k (λ) A p 1 (λ) = λ α 1 P = P 1 P 2 λ α 1 β 1 p 2 (λ) = β 1 λ α 2 P T AP = D AP = P D p k (λ) = λ α 1 β 1 0 β 1 λ α 2 β 2 0 0 β 2 λ α 3 β 3............... 0 0...... λ αk 1 β k 1 0 0 β k 1 λ α k (6.28) p k (λ) = (λ α k )p k 1 (λ) β 2 k 1p k 2 (λ) (6.29) k = 2 n n A p 0 (λ) = 1 p k(λ) = (λ α k )p k 1(λ) β 2 k 1p k 2(λ) +p k 1 (λ) (6.30) n p n (λ) p k (λ)

6. 8 3 p n (λ) (6.29) (6.30) p n (λ) n n A λ 1 > λ 2 > > λ n λ p n (λ), p n 1 (λ),, p 0 (λ) u k N(λ) N(λ) λ Au k = λ k u k (6.34) a < λ < b p n (a) 0 p n (b) (a, b) p n (λ) x 0 = c 1 u 1 + c 2 u 2 + + c n u n N(a) N(b) (6.31) r < λ < r lim ( ) (a, b) (6.33) c 1 0 x 0 u 1 x j A( ) λ 2 u 2 x 0 u 1 x j = Ax j 1 (6.32) j x j A λ 1 x 0 λ 1 lim j x T j x j x T j x = λ 1 (6.33) j 1 c 1 0 (6.32) x j = A j x 0 = c 1 λ j 1 u 1 + c 2 λ j 2 u 2 + + c n λ j nu n [ = c 1 λ j 1 u 1 + c ( ) j 2 λ2 u 2 + c 1 λ 1 (6.19) + c ( ) ] j n λn u n c 1 λ 1 r = max ( β i 1 + α i + β i ) 1 i n j x j x j = u 1 u 1

6. 9 A λ k µ µi A B = (µi A) 1 (6.35) λ k A u k B 1 u k = (µi A)u k = (µ λ k )u k 1 Bu k = u k µ λ k B (µ λ k ) 1 µ λ k λ i µ λ k < µ λ i i k B A QR (µ λ k ) 1 x 0 B A 1 = A = Q 1 R 1 A 2 = R 1 Q 1 lim j lim j x j x j = u k u k x T j x j x T j x = j 1 1 µ λ k (6.36) (µi A) 1 (µi A)x j = x j 1 (6.37) ( ) 2 2 x j 1 x j (6.36) j RQ A Q A (6.37) λ k. = µ x T j x j 1 x T j x j λ k µ QR (6.37) µ QR (6.37) LU QR µ µi A A λ v (6.37) A vi λ v A vi QR λ v λ v x j 1 = 1 µ λ k QR A QR A Q Q 1 AQ = Q 1 (QR)Q = RQ A RQ A 2 QR A k = Q k R k A k+1 = R k Q k (6.38) A k A A

6. 10 (6.40) A k v k I = Q k R k B [q 1, q 2, ] A k+1 = R k Q k + v k I (6.39) h 11 h 12 A k QR A k v k I QR = [q 1, q 2, ] h 21 h 22 0 h 32 (6.41) R k Q k = Q T k A k Q k v k I A k+1 = Q T k A k Q k Bq 1 = q 1 h 11 + q 2 h 21 A k+1 A k k A k q j QR v k q T 1 Bq 1 = h 11 v k A k q 1 (n, n) h 11 h 11 v k q 2 h 21 = Bq 1 q 1 h 11 n (n, n) λ n λ n n n h 2 21 = Bq (n 1) (n 1) 1 q 1 h 11 2 h 21 = Bq 1 q 1 h 11 v k h 21 q 2 = 1 (Bq h 1 q 1 h 11 ) 21 q 2 (6.41) B Bq 2 = q 1 h 12 + q 2 h 22 + q 3 h 32 BQ = QH (6.40) Q H h 12 = q T 1 Bq 2 h 22 = q T 2 Bq 2 Q Q H H h 12, h 22 Q H h 32 = Bq 2 q 1 h 12 q 2 h 22 Q q 3 = 1 (Bq h 2 q 1 h 12 q 2 h 22 ) 32 H Q h 32 q 3 Q q j H Q = [q 1, q 2,, q n ] q 1, q 2

6. 11 QR (6.9) (6.10) P k Q = P n 1 P n 2 P 1 A A Q QR P 1 A k v k I = Q T k R k A k+1 = R k Q k + v k I (6.42) A k+1 v k+1 I = Q T k+1r k+1 A k+2 = R k+1 Q k+1 + v k+1 I (6.43) QR Q T m 11 = p 1 = a 2 11 (v k + v k+1 )a 11 + a 12 a 21 R +v k v k+1 (6.42) (6.43) m 21 = q 1 = a 21 (a 11 + a 22 v k v k+1 ) (6.45) P 1 P 1 M 0 A k M A k a ij m 31 = r 1 = a 21 a 32 A k+1 = Q k A k Q T k A k+2 = Q k+1 A k+1 Q T k+1 A k+1 (6.43) s = sign( p 1 ) p 2 1 + q2 1 + r2 1 (6.46) A k v k+1 I = Q T k Q T k+1r k+1 R k (6.42) w 0 M (A k v k+1 I)(A k v k I) P 1 = Q T R (6.44) Q = Q k+1 Q k R = R k+1 R k Q A k+2 (6.10) w = 1 s(s p1 ) [p 1 s, q 1, r 1, 0,, 0] T A k+1 P 1 A k+1 A k P 1 A k A k+2 = QA k Q T A k+2 A k P 1 A k = 0 0 A k Q T = Q T 0 0 0 A k+2.... P 1 Q T 0 Q (6.44) P 1 M 0 Q M P 1

6. 12 p 2 q 2 P 1 A k P 1 = r 2 0 0 0........ q 2, r 2 0 P 2 1 0 n l a l,l 1 P 2 = 0 I n 1 ww T l = n a nn n n ww T (6.46) p 1, q 1, r 1 p 2, q 2, r 2 0 l n n 1 A k+2 (6.45) Q = P n 1 P n 2 P 1 Q a 11 a ll a 12 a l,l+1 a 32 a l+2,l+1 QR P k v k, v k+1 A k 2 2 a l,l 1 0 l λ a n 1,n 1 a n 1,n a n,n 1 λ a n,n = 0 (6.45) QR v k + v k+1 = a n 1,n 1 + a n,n v k v k+1 = a n 1,n 1 a n,n a n 1,n a n,n 1 A A k a n,n 1 a n,n A k (n 1) (n 1) n m (n m) A a n 1,n 2 A k l a l,l 1 0 l n l l 1 n l + 1 l = n 1 2 2 (1,1) (l, l) 2 2 A = UΛV T (6.47)

6. 13 U n m (6.47) U V V m m Λ m m U V U T U = I m V T V = I m B S T V V V T = I m U UU T U (6.47) U V (6.47) A Λ A U = P S V = QT 0 0 p A ( ) n m B p m p = m p < m A κ(a) = max(λ i) min(λ i ) (6.48) U T AV = Λ S T BT = Λ T 1 = 0 B A 0 0 r 3 0 B = BT = 0 0.. n m(n m) A 0 0 0 (6.9) P k 0 (6.17) 0 (2,1) 0 B Q k T 1 B(6.21) ϕ 1 θ 1 S 1 B = P T AQ (6.49) 0 P = P 1 P 2 P m B = S T B = 0 Q = Q 1 Q 2 Q m 1 m = 4 B w 1 r 2 0 0 B (2,1) 0 θ 1 0 w 2 r 3 0 (1,3) 0 0 0 w 3 r 3 B = (6.50) T 1 0 0 0 w 4... T 2 ϕ 2... 0 (1,3) 0 (3,2) 0 0 S T 2 0 cos ϕ 1 sin ϕ 1 0 sin ϕ 1 cos ϕ 1 0 0 I

6. 14 l w l 1 0 r l 0 B 2 1 4 3 6 5 T k S k I k 1 c s s c 7 I m k 1 l l = m S T m 1S T m 2 S T 1 BT 1 T 2 T m 1 T 1 ϕ 1 1.260 0.840 0.630 0.504 0.840 0.630 0.504 0.420 A = 0.630 0.504 0.420 0.360 T 1 0.504 0.420 0.360 0.315 [cos ϕ 1, sin ϕ 1, 0,, 0] T 0.420 0.360 0.315 0.280 B T B si 2.558200e+0 s 1.799580e 1 6.208599e 3 [w1 2 s, w 1 r 2, 0,, 0] T 9.967085e 5 s B T B 2 2 [ w 2 m 1 + r 2 m 1 w m 1 r m ] w m 1 r m wm 2 + rm 2 6 (ill-conditioned) QR B U V QR S T U S B T V T r m, w m 1, r m 1, w m 2, U V r l 0 w l 1 = 0 l l r l l 1 l r l (l 1, l) 0 (l 1, l + 1) l 1 l + 1 (l 1, l + 1) 0 (l 1, l + 2) l 1 k (l 1, k) l (6.48) κ(a) = 2.567 10 4 w m

7. 1 7 f(x) I = b a f(x)dx (7.1) h = b a ( a + b ) I C = hf 2 I I C = h2 4! [f (b) f (a)] (7.2) 7h4 2 3 6! [f (3) (b) f (3) (a)] f(x)dx = h + 31h6 x 0 2 (f 0 + f 1 ) 2 3 3 8! [f (5) (b) f (5) (a)] + h3 12 f (ξ) (7.5) x2 f(x)dx = h x 0 3 (f 0 + 4f 1 + f 2 ) h5 90 f (4) (ξ) (7.6) x3 I T = h f(x)dx = 3h x [f(a) + f(b)] (7.3) 0 8 (f 0 + 3f 1 + 3f 2 + f 3 ) 2 I I T = h2 3h5 2 3! [f (b) f (a)] 80 f (4) 3 (ξ) 8 (7.7) x4 + h4 6! [f (3) (b) f (3) (a)] h6 3! 7! [f (5) (b) f (5) (a)] + f(x)dx = 4h x 0 14175 (989f 0 + 5888f 1 928f 2 + 10496f 3 4540f 4 +10496f 5 928f 6 + 5888f 7 +989f 8 ) 2368h11 467775 f (10) (ξ) (7.9) h 2 f (b) f (a) hf (ξ) I S = h [ ( a + b ) ] f(a) + 4f + f(b) (7.4) 6 2 (7.4) (7.4) = 1 3 (2I b a = h (7.6) C + I T ) x 2 x 0 = 2h h4 I I S = 2 2 6! [f (3) (b) f (3) (a)] + 5h6 2 3! 8! [f (5) (b) f (5) (a)] + f(x) x f j = f(x j ) h = x j+1 x j ξ x1 x 0 x8 f(x)dx = 2h 45 (7f 0 + 32f 1 + 12f 2 +32f 3 + 7f 4 ) 8h7 945 f (6) (ξ) (7.8) 2

7. 2 h x 0 x n x 0 x n f(x)dx = 3h f(x) 2 (f 1 + f 2 ) + h3 4 f (ξ) f(x) = e iωx 4 (7.8) h(iωh)6 118.1 eiωξ f(x) 8 (7.9) (7.3) (7.4) h(iωh)10 197.5 eiωξ ωh > 1 (7.9) cos 2 θdθ = π 0 2 f(x) h = π/2 θ = 0, π/2, π 10 h ωh < 2π/10 < 1 I T = π ( 1 2 2 + 0 + 1 = 2) π 2 x2n f(x)dx = h [ 1 x 0 2 f 0 + f 1 + f 2 + +f 2n 1 + 1 2 f ] 2n (7.10) h2 [ f (x 2n ) f (x 0 ) ] + 12 (7.4) h 4 h 6 x2n f(x)dx = h x 0 3 [f h 0 + 4f 1 + 2f 2 + +2f 2n 2 + 4f 2n 1 + f 2n ] (7.11) h4 180 [ f (3) (x 2n ) f (3) (x 0 ) ] + x3 x 0 f (2k 1) (x) 0 π (7.9) h (7.8) h h x 2n 1 1 1 + x 2 dx 1 f (x) = 24x(1 x2 ) (1 + x 2 ) 4 [a, b] n h = b a n x j = a + jh j = 0, 1, 2,, n

7. 3 h I T (h) = h [ 1 2 f 0 + f 1 + + f n 1 + 1 2 f ] n (7.12) h h/2 h/4 h = 0 m = 5 I I T (h) = αh 2 + βh 4 + α, β h I I T (h/2) = 1 4 αh2 + 1 16 βh4 + h 2 (h = 0.5) I = 1 3 [4I T (h/2) I T (h)] 1 4 βh4 + k S (k) 0 S (k) 1 S (k) 2 S (k) 3 (7.9) (7.10) 0 9.254511 8.155954 8.153362 8.153365 1 8.430593 8.153524 8.153365 h 4 2 8.222792 8.153375 3 8.170728 S (k) S (0) 0 = I T (h/2 k 3 ) k = 0, 1, 2, 32 S (k) m S (k) m = S (k+1) m 1 + S(k+1) m 1 S(k) m 1 4 m 1 m = 1, 2, (7.13) m k e x2 dx S (0) 0 S (0) 1 S (0) 2 S (1) 0 S (1) 1 S (2) 0 m 2 S (0) 0 h S (1) 0 S (0) h = 1 1 h S (2) 0 S (1) 0 S (2) 0 S (1) x = 0, ±1, ±2, ±3, ±4 1 S (0) 1 S (0) 1.772637 π = 1.772453 2 h = 0.5 10 S m (k) S (1) 0 S (0) 0 S (0) 1 S (0) 0 S (0) 1 1 I = f(x)dx S (1) 1 S (1) 0 1 S (0) 2 S(1) 1 S(2) 0 x h x = tanh y dx = dy f(x) cosh 2 (7.14) y 0 x 4 log(x + 1 + x 2 )dx (7.3) 0 f(x) = e x2 x = ± 0

7. 4 y I = dy f(tanh y) cosh 2 y y f(tanh y)/ cosh 2 y y = ± e 2 y 0 π dx = 2 cosh y ( cosh y 2 π )dy 2 sinh y h I SE = h w k f(x k ) (7.15) y ± k= 1 y k = kh x k = tanh y k w k = cosh 2 y k ± w k 0 (DE) h 1/2 1/4 h I DE (h) = h w k f(x k ) (7.18) x k w k k= ( h = 1/2 M π ) y k = kh x k = tanh 2 sinh y k π h = 1/2, 1/4, 1/8, x k w k w k = 2 cosh y ( k (7.15) h cosh 2 π ) 2 sinh y k w k (7.16) h I SE (h) ( I I SE (h) exp c ) (7.16) h c f(t) I I SE (h/2) (I I SE (h)) 2 I SE (h/2) I I I SE (h/2) I ISE (h/2) I SE (h/2) I SE (h) h t f(x) = 1 2 (7.14) ( π ) x = tanh 2 sinh y (7.17) ( 1975) y exp( c exp y ) F (p) = 0 f(t)e pt dt p > 0 p e pt t p pt = e y dt = 1 p ey dy (7.19) h k w k = 2 F (p) = 1 f(t)w(y)dy (7.20) p

7. 5 w(y) = exp(y e y ) (7.21) n 2 ±x k 0.57735 02692 n 3 ±x k 0 0.7071067812 4 0.18759 24741 5 0 F (p) = h f(t k )w k (7.22) 0.79465 44723 0.37454 14096 p k= 0.83249 74870 y k = kh pt k = e y k w k = w(y k 6 0.26663 54015 7 0 0.42251 86538 0.32391 18105 (7.16) 0.86624 68181 (7.16) 0.52965 67753 0.88386 17008 y = 0 (7.19) 9 0 w(y) 0.16790 61842 10 10 0.52876 17831 y = 20.7 pt = 1.0 10 9 0.60101 86554 0.91158 93077 y = 2.0 pt = 20.7 e pt pt > 20.7 t = 0 1 [ 1, 1] 1 n I = f(x)dx I G = w k f(x k ) (7.24) 1 I Ch = 2 n f(x k ) (7.23) n k=1 f(x) n x 2n 1 I Ch 0 x k n = 2 n = 3 k=1 2n x k, w k f(x) x 1 + x 2 = 0 x 2 1 + x 2 2 = 2 3 w 1 + w 2 + w 3 = 2 x 2 = x 1 = 1 w 1 x 1 + w 2 x 2 + w 3 x 3 = 0 3 w 1 x 2 1 + w 2 x 2 2 + w 3 x 2 3 = 2 3 w 1 x 3 1 + w 2 x 3 2 + w 3 x 3 3 = 0 n = 8 n > 9 w 1 x 4 1 + w 2 x 4 2 + w 3 x 4 3 = 2 5 w 1 x 5 1 + w 2 x 5 2 + w 3 x 5 3 = 0

7. 6 x 2 = 0 x 1 = x 1 x 2 1 = 3 5 n ±x k w k w 1 = w 3 = 5 w 2 = 8 2 0.57735 02692 1 9 9 3 0 0.88888 88889 0.77459 66692 0.55555 55556 4 0.33998 10435 0.65214 51549 0.86113 63116 0.34785 48451 5 0 0.56888 88889 w(x) [a, b] p j (x) 5 p n (x) x k (k = 1, 2,, n) n 1 [ n 1 n ] 1 f n 1 (x) = w k p j (x k )f(x k ) p j (x) λ j=0 j k=1 f(x) f n 1 (x) I G = b a n 1 = λ j=0 j w(x)f n 1 (x)dx [ n ] 1 b w k p j (x k )f(x k ) w(x)p j (x)dx k=1 b a w(x)p j (x)dx = 1 µ 0 b a w(x)p 0 (x)p j (x)dx = λ 0 µ 0 δ 0j j j = 0 I = b n a w(x)f(x)dx. = I G = k=1 a w k f(x k ) (7.25) w(x) = 1 x = y 2 p n (x) P n (x) P n (x k ) = 0 k = 1, 2,, n w k = 2(1 x2 k ) [np n 1 (x k )] 2 (7.26) n 1 x = 0 2n 1 0 0.53846 93101 0.47862 86705 0.90617 98459 0.23692 68851 f(x)e x dx. = k w k f(x k ) 1 g(x) dx x 0 g(x) x = 0 I = 2g(y 2 )dy 0

7. 7 I = g(x) xdx 0 x = 0 x x = 0 x = 0 (7.5) (7.6) (7.8) (7.9) I = 2g(y 2 )y 2 dy 0 h I(h) h/2 I(h/2) I = g(x)x ±α dx 0 < α < 1 h 0 I(h) 4 (7.8) x = 0 I I(h) = 8h7 945 f (6) (ξ) (7.16) y = ± h/2 I(h/2) x x = ±1 I I(h/2) = 8h7 2 7 945 [f (6) (ξ L ) + f (6) (ξ R )] ξ L ξ R x = 1 f (6) (ξ) (7.16) exp ( π f ) (6) (ξ) 1 x = 2 sinh y ( π ) I I(h/2) = 1 [I(h/2) I(h)] (7.27) cosh 2 sinh y 63 (7.17) f(x) ε x 1 x 1 I(h/2) I(h) < ε (7.28) x = 1 f(x) 63 x = 1 I(h/2) I(h/2) (7.23) I =. I(h/2) + 1 [I(h/2) I(h)] (7.29) 63 I(h/2) I(h/4)

7. 8 ε ε/2 0 1 2 1 2 3 1 4 0 1 2 3 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 5 8 9 9 10 5 5 10 1975

8. 1 8 y n = Az1 n + Bz2 n (8.4) A B n = 0 y 0 y 1 y 0 = A + B y 1 = Az 1 + Bz 2 (8.5) d 2 y(t) dt 2 + ω 2 y(t) = 0 (8.1) t = n t y n = y(n t) z 1 d 2 y(t) dt 2. = y n+1 2y n + y n 1 t 2 (8.3) A B z 1 < z 2 B 0 n z 1 A = 0 N y N y N+1 A B (8.3) n = N, N 1,, y n+1 [2 (ω t) 2 ]y n + y n 1 = 0 (8.2) y n 1 = 1 β (αy n + y n+1 ) n y n+1 + αy n + βy n 1 = 0 (8.3) z 2 (8.3) α β n z 1 < z 2 (8.3) y 0 y 1 z 1 n = 0, 1, y n z 2 z 1 = z 2 α β n y n z 1 < 1 (8.3) n y n y n z n (8.3) n + z n z n (8.3) z 1 > 1 z 2 > 1 (8.6) (8.3) z 2 + αz 1 + β = 0 z z z1 1, z 1 2 = 1 [ α ± ] ω 2 > 0 (8.2) α 2 2 4β 1 z 1

8. 2 (8.2) ω t 1.. z 1, z 2 = 1 ± iω t = e ±iω t d 2 Z n (x) dx (8.2) (8.1) 2 + 1 ) dz n (x) + (1 n2 x dx x 2 Z n (x) = 0 (8.8) y(t) = Ae iωt + Be iωt (8.3) α = 2 cos θ β = 1 z 1 1, z 1 2 = cos θ ± i sin θ = e ±iθ 1 (8.3) x J n (x) y 0 = 0 x A + B = 0 J 0 (x) x = 10 y 1 = A(e iθ e iθ ) = 2iA sin θ k a k S k A = 1 0 1.000000e+0 1.000000e+0 2i 1 2.500000e+1 2.400000e+1 (8.3) 2 1.562500e+2 1.322500e+2 y n = sin nθ y 0 = 1 y n = cos nθ 8 9.385967e+1 2.275143e+1 9 2.896903e+1 6.217607e+0 (8.3) 10 7.242259e+0 1.024651e+0 sin(n + 1)θ + sin(n 1)θ = 2 cos θ sin nθ cos(n + 1)θ + cos(n 1)θ = 2 cos θ cos nθ (8.7) sin θ J 0 (10) cos θ cos nθ sin nθ 0.2459 3 J 0 (x) Z n (x) n 0 x 0 ( x ) n ( 1) k ( x ) 2k J n (x) = (8.9) 2 k!(k + n)! 2 k=0 3 4.340278e+2 3.017778e+2 4 6.781684e+2 3.763906e+2 5 6.781684e+2 3.017778e+2 6 4.709503e+2 1.691725e+2 7 2.402808e+2 7.110825e+1 a k (8.9) S k k k = 4 376.4 x = 5

8. 3 J n (x) Z n (x) Z n+1 (x) 2n x Z n(x) + Z n 1 (x) = 0 (8.10) Z 0 (x) Z 1 (x) y Z 2 (x), Z 2 (x), N = ε y N+1 = 0 (8.13) ε Z n (x) J n (x) (8.10) J 0 (1) J 1 (1) n J n (1) 0 0.765197686 1 0.440050585 2 0.114903484 3 0.019563351 4 0.002476622 5 0.000249625 y n J n (x) n = 5 3 (8.10) (8.3) α n J k (x) J k (x) =. y k (8.15) z1 1, z 1 2 = n S n x ± 2 x 2 1 ε n x n 1 (8.14) (8.10) n y k y k S n x (8.10) 1 (8.10) N n (x) n (8.10) Z n (x) z1 n z 1 < 1 N 0 (x) N 1 (x) J n (x) J 0 (x) + 2 J 2 0 (x) + 2 J 2n (x) = 1 (8.11) n=1 Jn(x) 2 = 1 (8.12) n=1 J n (x) N n N y k 1 = 2k x y k y k+1 (8.14) S = y 0 + 2 k = N, N 1,, 1 N 1 (x) n J n (x) 0 (8.10) N k=1 y 2k (8.11) S J 0 (x) J 1 (x) (8.13, 14) J n (x)n n+1 (x) J n+1 (x)n n (x) = 2 πx (8.16) J 0 (x) J 1 (x) N 0 (x)

8. 4 x < 5 N 0 (x) K n (x) π 2 N 0(x) = (ln x ) 2 + γ K 0 (x) K 1 (x) J 0 (x) K n (x) + 1 ( x ) 2 1 ( x ) 4 (1 + 1/2) (1!) 2 2 (2!) 2 2 1 ( x ) 6 +(1 + 1/2 + 1/3) I n (x)k n+1 (x) + I n+1 (x)k n (x) = 1 (8.22) (3!) 2 2 x K 0 (x) γ K n (x) x γ = 0.5772156 K 0 (x) N 0 (x) d 2 Z n (x) dx 2 + 1 dz n (x) x dx (1 + n2 x 2 ) Z n (x) = 0 (8.17) j 0 (x) = x 1 sin x I n (x) j 1 (x) = x 2 (sin x x cos x) (8.23) I n (x) = i n J n (ix) (8.18) I n (x) (8.9) I n (x) I n+1 (x) + 2n x I n(x) I n 1 (x) = 0 (8.19) j n (x) j n 1 (x) + j n+1 (x) = 2n + 1 j n (x) x (8.24) (8.23) z1 1, z 1 2 = n n n x ± 2 x 2 + 1 J n (x) 1 J n (x) (8.24) n I n (x) (8.19) j 0 (x) (8.23) I 0 (x) + 2 I n (x) = e x (8.20) n=1 j n (x) x j 2 = x 3 [(3 x 2 ) sin x 3x cos x] j 0 (x) x = kπ 0 j 1 (x) n n (x) j n (x) (8.24) N n (x) (8.17) K n (x) K n+1 (x) 2n x K n 0 (x) = x 1 cos x n(x) K n 1 (x) = 0 (8.21) n 1 (x) = x 2 (cos x + x sin x) (8.25) I n (x) n 2 (x) = x 3 [(3 x 2 ) cos x + 3x sin x] n n

8. 5 Z 0 S c n c n c k = a k P 0 (x) = 1 P 1 (x) = x x 1 (n + 1)P n+1 (x) (2n + 1)xP n (x) +np n 1 (x) = 0 (8.26) ( P n+1 (x) 2 1 ) xp n (x) n + 1 ( + 1 1 ) P n 1 (x) = 0 n + 1 n Z k x = cos θ (8.7) (8.23) Z k = U k + V k sin θ (8.29) (c) c n S = N c n e inθ (8.27) n=0 V N = 0 a N z = e iθ V 0 sin θ S S z N U 0 (a) k = 0 S = c 0 + z(c 1 + z(c 2 + + z(c N ) ) n Z n (8.30) k = 0 S = Z 0 = c 0 + zz 1 Z 1 = c 1 + zz 2 Z k = c k + zz k+1 a k (8.30) V 0 V 1 Z N = c N N N a k cos nθ a n sin nθ Z N+1 = 0 Z k = c k + zz k+1 k = N, N 1,, 0 (8.28) (8.28) Z k+1 = z 1 (Z k a k ) Z k 1 = zz k + a k 1 Z k 1 + Z k+1 = (z + z 1 )Z k + a k 1 z 1 a k V N = 0 V N+1 = 0 (a) (b) (c) V k 1 + V k+1 = 2 cos θv k + a k (8.30) k = N, N 1,, 1 U 0 = V 1 + V 0 cos θ + a 0 (8.31) U 0 = V 1 V 0 cos θ n=0 n=0

8. 6 p n+1 (x) = α n p n (x) + β n 1 p n 1 (x) (8.32) p n (x) S = N a n p n (x) (8.33) n=0 V 1, V 2, V 3 (8.32) V 0, V 1 V 2 p n+1 (x) a n+1 sin nθ cos nθ (8.7) Fig. 8.1 1 a n V 0 sin θ 1 1 a k n S n = a k V k V n = sin nθ (d) k=0 V n V n 1 + V n+1 = 2 cos θv n (e) S n sin cos V n+1 (8.31) S n+1 = S n + a n+1 V n+1 (f) Fig. 8.1 (8.32) (8.33) Fig. 8.2 q n p n (x) α p n 2 cos θ β 1 sin, cos β β q N+1 = q N+2 = 0 1 V 1 V 2 sin θ α V 3 α q n = a n + α n q n+1 + β n q n+2 (8.34) n = N, N 1,, 1 a 1 a 2 a 3 q 1 q 0 = a 0 + β 0 q 2 S 1 S 2 Fig. 8.1 S 3 S N (8.30) Fig. 8.1 (d) sin nθ cos nθ Fig. 8.1 V 1 V 1 S = q 0 p 0 (x) + q 1 p 1 (x) (8.35) 1 V 1, V 2 a n

8. 7 1 p 1 (x) p 0 (x) q 0 β 0 β 1 q 1 q 2 α 1 β N 2 q N 1 αn 1 q N a 0 a 1 a 2 a N 1 a N S 0 S 1 S 2 S N 1 S N Fig. 8.2 x n (8.3) x n = x 0 δ n x 0 y n + αy n 1 + βy n 2 = x n (8.36) y y n = x 0 h n n = 0 n = k x k x n = 0 y n = 0 n < 0 k (8.36) y n = x k h n k n 0, 1, 2, y n h n x 0, x 1, x 2 { n n 1 n = 0 h n +αh n 1 +βh n 2 = δ n = 0 n 0 (8.37) y n = x k h n k = h k x n k (8.39) k=0 k=0 δ n h n x n h n h n n 0 n y n h k x n k max x k h k k=0 k=0 { 1 n = 0 (8.38) h n = 0 n < 0 n (8.6) (8.6) (8.36) y n = x n αy n 1 βy n 2 (8.40) h n n = 0, 1, 2, h n = αh n 1 βh n 2 n = 1, 2, (8.39) h n < (8.38) (8.39) (8.40) n=0 (8.4) (8.6)

9. 1 9 t x j x j 2 < (9.1) j= x j ( 14) x j X D (ω) = x j e ijω (9.2) j= x j ω f x j t ω = 2πf t (9.3) X D (ω) X D (ω + 2π) = X D (ω) (9.4) ω π (9.5) f f N = 1 (9.6) 2 t 2 t t 2 t x j (9.1) x j < (9.2) (9.1) π N x j e ijω X D (ω) 2 dω = 0 lim N π j= N X D (ω) (9.2) x j = 1 π X D (ω)e ijω dω (9.7) 2π π (9.2) x j e ijω (9.3) x j = fn f N tx D (ω)e ijω df tx D (ω) (9.2) (9.7) 5 (5.16) x f(x) (a k, b k ) x j x(t) t t x(t) (aliasing) Fig. 9.1 ( ) x(t) X(σ) = x(t)e iσt dt (9.8)

9. 2 σ σ = 2πf x(t) = 1 2π X(σ)e iσt dσ (9.9) X D (ω) X(ω) x j = x(j t) (9.2) j= x(j t)e ijω x(t) (9.9) f(t) F (σ) αβ = 2π α, β β f(mα) = F (nβ) (9.10) 2πα m= n= f(t) f(t) = x(t)e iωt/ t j f(j t) f(t) F (σ) = X(σ + ω/ t) α = t X D (ω) = = 1 t m= n= x(m t)e imω ( ) ω + 2πn X t X(σ) X D (ω) ( ) ω + 2πn tx D (ω) = X (9.11) t n= X D (ω) ω X(σ) σ X Fig. 9.2 X(ω) X D (ω) X(ω) ω > π π < ω X(ω) ω < π π < ω X(ω) f X D (ω) ω > π X( ) 0 X D ( ) 0 Fig. 9.2 f f f, f ± 1/ t, f ± 2/ t,, f ± n/ t, f ± 2nf N f Fig. 9.1 f = 13 Hz t = 0.1 s f N = 5 Hz f = 13 2 5 = 3 Hz 3 Hz X(σ) 0 X D (ω) X(σ) tx D (ω) = X(ω/ t) (9.12)

9. 3 σ > π X(σ) = 0 t X D (ω) D X(ω) x j, y j z j = x j y j (9.13) Z(ω) Z(ω) = 1 2π π π X(σ)Y (ω σ)dσ (9.14) z j x j y j X Y z j = x k y j k = x j k y k (9.15) k= k= z j Z(ω) = X(ω)Y (ω) (9.16) y j = x j (9.15) { } z j = k x k x k j = k x j+k x k x k y j Y (ω) = X(ω) z j Z(ω) = X(ω) 2 z j = 1 2π π π X(ω) 2 e ijω dω j = 0 z 0 = x k 2 = 1 π X(ω) 2 dω (9.17) 2π k= π x j t X(ω) 2 (9.2) (9.7) (9.8) (9.9) N x 0, x 1,, x N 1 c k c k = N 1 j=0 ( ) 2πijk x j exp N k = 0, 1, 2,, N 1 x j = 1 N N 1 k=0 ( c k exp 2πijk ) N (9.18) (9.19) (9.18) N 1 [ ] { 2πik(l j) N l = j exp = N 0 l j k=0 (9.18) (9.2) ω k = 2π N k c k = X(ω k ) (9.20) X(ω) ω = 2π N

9. 4 c k (9.7) (9.19) X(ω) (9.7) (0, 2π) N 1 2π 2π 0 X(ω)e ijω dω. = ω 2π = ω 2π N k=0 N 1 k=0 X(ω k )e ijω k X(ω k )e ijω k k = 0, N 1/2 X(ω) X(0) = X(ω N ) X(ω k ) = c k (9.19) (9.7) ω N N N N = 2M + 1 2M + 1 x j ( j M) M X(ω) = x j e ijω j= M c k = X(ω k ) ω k = 2πk 2M x j = 1 2M M k= M k M (9.21) c k e ijω k (9.22) j = M j = M 1/2 (9.21) c k X(ω) (9.7) ( π, π) 2M ω (9.15) x j, y j j = 0 N 1 N z j (9.15) (9.16) z j = 1 2π 2π 0 X(ω)Y (ω)e ijω dω FFT x j, y j N z j ~ z x z Fig. 9.3 x, y N z x y x y N z j j = 0 2N 2 N z j z j N z j+n = z j z j j = N 2N 1 j = 0 N 1 z j = n= z j+nn z j z j 2N y

9. 5 x j (9.18) (9.19) N N N x j (9.18) c k c N k = c k N c k k = 0 N/2 c k c k = a k + ib k a k, b k a k = b k = N 1 j=0 N 1 j=0 x j cos jω k ω k = 2π N k x j sin jω k (9.23) k = 0, 1, 2,, N/2 b 0 b N/2 0 a k, b k N N (9.19) k N 1 k=0 c k = N 1 k=n/2 N/2 1 k=0 N/2 c k = = N/2 k=1 k=1 c k exp c k + [ c N k exp N 1 k=n/2 ( 2πijk ) N c k ] 2πij(N k) N x j = 2 [ N/2 1 1 N 2 a 0 + (a k cos jω k + b k sin jω k ) k=0 + 1 ] 2 a N/2 cos 2πj (9.24) N k N/2 a N/2 N (0, 2π) N 2M +1 (9.21) a k = b k = M j= M M j= M x j x j cos jω k ω k = 2πk 2M x j sin jω k (9.25) x j = 1 [ M 1 1 M 2 a 0 + (a k cos jω k + b k sin jω k ) k=1 + 1 ] 2 a M cos jπ (9.26) j = + x j (9.2) X M (ω) = M j= M x j e ijω X(ω) X M (ω) x j (9.7) π X M (ω) = 1 M X(σ) e ij(ω σ) dσ 2π π j= M = 1 π X(σ)W M (ω σ)dσ (9.27) 2π π = 1 π X(ω σ)w M (σ)dσ 2π π W M (ω) = M j= M e ijω = sin(n + 1/2)ω sin ω/2 (9.28) (9.27) W M (ω) Fig. 9.4 X M (ω)

9. 6 W M (σ) ω X(ω σ) W M (σ) X M (ω) X(ω) W M (σ) W M (σ) ( ) ( ) 20 X(ω) 5 (Bartlett) w j = 1 j M W M (ω) = M j M ( ) 2 sin Mω/2 (9.31) M sin ω/2 W M (ω) 0 ω = 2π/M (9.31) W M (ω) Fig. 9.5 0 Box-Car M = 10 Fig. 9.4 j > M 0 w j { w j x j j M z j = (9.29) 0 j > M (9.14) z j Z(ω) = 1 2π π π W M (σ)x(ω σ)dσ (9.30) W M (ω) w j X(ω) z j x j Z(ω) X(ω) w j W M (ω) { 1 j M w j = 0 j > M ( j > M w j = 0 ) Cosine-bell ( ) ( ) w j = 1 ( 1 + cos πj ) j M 2 M sin Mω cos ω/2 W M (ω) = 2 sin ω/2 [ ( ) ] 2 1 sin ω/2 1 (9.32) sin π/2m Mω = π ω = 2π/M (Fig. 9.5 ) a + b ( 1 + cos πj ) j < M w j = 2 M 1 2 a j = M sin Mω cos ω/2 W M (ω) = sin ω/2 { [ a + b ( ) ] 2 1 } sin ω/2 1 (9.33) 2 sin π/2m a, b ω = 5π/2M

9. 7 0 a, b W M (5π/2M) = 0 w 0 = 1 a, b M > 10 M a = 4 46 b 2 = 21 46 (9.34) ( ) (Fig. 9.5 ) M = 10 0 Triangle (Bartlett) M = 10 (9.28) (Fig. 9.4 ) 1 j < M w j = 1 j = M 2 W M (ω) = sin Mω cos ω/2 sin ω/2 (9.35) j 0 W M (ω j ) = 0 j = 0 ( ) 2 3 j 1 6 + 6 j M M w j = 0 j M/2 ( 1 j ) 3 M/2 j M M 0 Cos-bell (Hann) M = 10 0 Hamming W M (ω) = 3 ( ) 4 sin Mω/4 4 M (9.36) M/2 sin ω/2 w j W M (ω) M (9.31) Fig. 9.5 Fig. 9.5 X(ω) X(ω) (9.30) Z(ω k ) = 1 2M M j= M W M (ω j )X(ω k j ) (3.31) (3.32) (3.33) W M (ω j ) j > 2 0 j = 0, ±1 FFT x j k a k, b k 2N k N N 2 N 3 N N N 2 1965 Cooley and Tukey (FFT, Fast Fourier Transform) FFT

9. 8 2N log 2 N N = 2 10 = 1024 N/2 log 2 N = 1024/20 = 50 FFT N N = L M j k j = l + ml k = n + pm 0 l < L, 0 m < M (9.37) n < M, 0 p < L (l, m) j x j L M x j x(l, m) c k M L c(n, p) x(l, m) c(n, p) jk N = ln N + lp L + mn M + mp e(x) = exp(2πix) (9.18) ( jk ) ( ln ) ( lp ) ( mn ) e = e e e N N L M mp e(mp) = 1 (9.18) L 1 ( lp )( ln ) c(n, p) = e L N l=0 M 1 ( mn ) e x(l, m) (9.38) M m=0 M 1 m=0 ( mn e M ) x(l, m) l M x(l, m) l 0 n < M M c (l, n) = M 1 m=0 ( mn e M ) x(l, m) c (l, n) L M l M x l c (l, n) l c (l, n) (9.38) L 1 ( lp ) ( ln ) c (p, n) = e e c (l, n) L N l=0 c (l, n) n e(ln/n) L c k = c(n, p) c(n, p) = c (p, n) c c FFT c (l, m) M 2 c (p, n) L 2 M 2 L + L 2 M + N = N(L + M + 1) (9.18) N 2 L M N N = n 1 n 2 n m N(n 1 + n 2 + + n m ) N = 2 m 2mN = 2N log 2 N (9.39)

9. 9 N 2 L = 2 ( n ) c (0, n) = c (0, n) + e c (1, n) ( N n ) c (1, n) = c (0, n) e c (1, n) (9.40) N c c (0, n) x j c (1, n) c (9.40) N = 2 3 j, k j j = (j 2 j 1 j 0 ) 2 = j 2 2 2 + j 1 2 1 + j 0 2 0 ( ) 2 k jk = j 0 (k 2 k 1 k 0 ) 2 + 2j 1 (k 1 k 0 ) 2 + 2 2 j 2 k 0 mod(n) mod(n) N N 1 ( ) jk e x(j) = ( j0 k ) ( 2 j0 (k 1 k 0 ) ) 2 e e N 2 N j=0 j 0 ( j1 k ) ( 1 2j1 k ) 0 e e 2 N j 1 ( j2 k ) 0 e x(j 2 j 1 j 0 ) (9.41) 2 j 2 j j 0, j 1, j 2 0 1 j 2 k 0 = 0 : x(0j 1 j 0 ) + x(1j 1 j 0 ) k 0 = 1 : x(0j 1 j 0 ) x(1j 1 j 0 ) x x(0j 1 j 0 ), x(1j 1 j 0 ) x(k 0 j 1 j 0 ) = ( j2 k ) 0 e x(j 2 j 1 j 0 ) 2 j 2 j 1 j 0 (j 1 j 0 ) 2 0 (11) 2 = 3 j 1 j 0 x(k 0 k 1 j 0 ) = ( j1 k ) ( 1 2j1 k ) 0 e e x(k 0 j 1 j 0 ) 2 N j 1 x(k 0 k 1 k 2 ) = ( j0 k ) 2 e 2 j 0 ( ) j0 (k 1 k 0 ) 2 e x(k 0 k 1 j 0 ) N x(k 0 k 1 k 2 ) c(k) c(1) = c(001) x x(100) = x(4) x(k 0 k 1 k 2 ) c(k 2 k 1 k 0 ) (k 1 k 0 ) 2 x x(k 0 k 1 j 0 ) (k 1 k 0 ) x (k 0 k 1 ) k 0, k 1 x(j) x(j 0 j 1 j 2 ) x(j 0 j 1 k 0 ) = ( j2 k ) 0 e x(j 0 j 1 j 2 ) 2 j 2 x(j 0 k 1 k 0 ) = ( j1 k ) ( 1 2j1 k 0 e e 2 N j 1 x(k 2 k 1 k 0 ) = ( j0 k ) 2 e 2 j 0 ( j0 (k 1 k 0 ) 2 e N ) x(j 0 j 1 k 0 ) ) x(j 0 k 1 k 0 ) (9.42) x(k 2 k 1 k 0 ) c(k) Cooley-Tukey Cooley et al. (1969)

9. 10 Subroutine FFT(A, M, NML) Complex A(1024), U, W, T Data PI/3.14159265/ SPI = PI If( NML.lt.0 ) SPI =-PI N = 2**M NV2 = N/2 J = 1 Do I=1, N-1 If( I.lt.J ) then T = A(J) A(J) = A(I) A(I) = T Endif K = NV2 6 If( K.ge.J ) Go to 7 J = J - K K = K/2 Go to 6 J = J + K Enddo Do L=1, M LE = 2**L LE1 = LE/2 U = 1 W = Cmplx(Cos(PI/LE1), + Sin(SPI/LE1)) Do J=1, LE1 Do I=J, N, LE IP = I + LE1 T = A(IP)*U A(IP) = A(I) - T A(I) = A(I) + T Enddo U = U*W Enddo Enddo Return End Cooley et al. (1969) FFT NML (9.18) 1 (9.19) 1 (9.19) 1/N Do Do j k jk = k 0 (j 2 j 1 j 0 ) 2 + 2k 1 (j 1 j 0 ) 2 + 2 2 k 0 j 0 mod(n) N 1 ( jk ) e x(j) = ( k2 j ) 0 ( k1 j ) 1 e e N 2 2 j=0 j 0 j 1 ( 2k1 j ) 0 ( k0 j ) 2 e e N 2 j 2 ( (j1 j 0 ) 2 k ) 0 e x(j 2 j 1 j 0 ) (9.43) N Sande-Tukey 0 N/2 b(k) = N 1 j=0 [x(2j) + ix(2j + 1)] e ( ) 2jk N k = 0, 1,, N/2 (9.44) c(k) = N/2 1 +e j=0 ( k N x(2j)e ) N/2 1 j=0 ( ) 2jk N x(2j + 1)e ( ) 2jk N x(j) c(k) = 1 {[ ] b(k) + b(n/2 k) 2

9. 11 ( ) k [ ie b(k) b(n/2 k)]} N k = 0, 1,, N/2 c(k) N c(n k) = c(k) (9.45) c(n/2 + k) = c(n/2 k) x(j) = 1 N ( +e j 2 N/2 1 k=0 [ c(k) ) c(n/2 k) ] e ( jk ) N b(k) = [ c(k) + c(n/2 k) ] ( +ie k ) [c(k) ] c(n/2 k) (9.46) N k = 0, 1,, N/2 x(2j) + ix(2j + 1) = 1 N N/2 1 e k=0 b(k) ( 2jk ) N (9.47) N/2 Cooley, J. W. and J. W. Tukey (1965) : An algorithm for the machine computation of complex Fourier series, Math. Comput., 19, 297-301. Cooley, J. W., P. A. W. Lewis, and P. D. Welch (1969) : The fast Fourier transform and its applications, IEEE Trans., Education, 12, 27-34.

10. 1 10 b a j b = (t j t )(y j y) j (t a = y bt (10.4) j t ) 2 b t j y j j r = (t j t )(y j y) j (t j t ) 2 j (y j y) 2 t j y j t j y j (t j, y j ) t 1 t y n y j a + bt j (10.1) a, b S = n (y j a bt j ) 2 = min (10.2) j=1 a b S a = 2 (y j a bt j ) = 0 S b = 2 (y j a bt j )t j = 0 a j a j 1 + b j t j + b j t j = j t 2 j = j y j t j y j (10.3) a, b t j y j t = 1 n n j=1 t j y = 1 n n j=1 y j (10.1) a b n y j y = [y 1 y 2 y n ] T T m x j x = [x 1 x 2 x m ] T y i x j y i y i a i1 x 1 + a i2 x 2 + + a im x m a ij n m a 11 a 12 a 1m a 21 a 22 a 2m A =...... a n1 a n2 a nm y = Ax + r (10.5) r

10. 2 x = [a, b] T 1 t 1 1 t 2 A =.. 1 t n r L 2 (10.6) S = r 2 = (y Ax) T (y Ax) (10.7) x δs = δx T A T (y Ax) (y Ax) T Aδx S δs δx δx 0 A T Ax = A T y (10.8) A T A m m A T y m m x = (A T A) 1 A T y (10.9) A T A A 1 1 1 A = ε 0 0 0 ε 0 (10.10) 0 0 ε A T A = 1 + ε 2 1 1 1 1 + ε 2 1 1 1 1 + ε 2 (A T A) 1 1 = ε 2 (3 + ε 2 ) 2 + ε 2 1 1 1 2 + ε 2 1 1 1 2 + ε 2 ε (A T A) 1 [1, 1, 1] T = 1 [1, 1, 1]T 3 + ε2 ε 2 1 (10.8) A 5 4 1.26 0.84 0.63 0.504 0.84 0.63 0.504 0.42 A = 0.63 0.504 0.42 0.36 0.504 0.42 0.36 0.315 0.42 0.36 0.315 0.28 y x x = [1/2, 1/3, 1/4, 1/5] T (10.11) y = Ax A T A A T y (10.8) 0.41049 x = 0.95418 0.92842 0.86204 A T A 10 6 damped least squares x S = min y Ax 2 + λ x 2 = min

10. 3 x λ r x λ [ ] [ ] y A = x + r (10.12) 0 λim r n + m (A T A + λi m )x = A T y A T A λ λ 1 + λ (A T A) ii (1 + λ)(a T A) ii (10.11) λ = 0.001 0.46900 x = 0.37860 0.27395 0.16619 λ (10.8) A T (y Ax) = A T r = 0 (10.13) S r A A a j Ax = a 1 x 1 + a 2 x 2 + + a m x m y a 1, a 2,, a m r r a j a T j r = 0 j = 1, 2,, m (10.13) a j a 1 y x 1 y a 1 a 1 x 1 a 2 a 1 y x 2 a 2 a j a j 6 ( ) A QR A = QR Q n m R m m Q Q T Q = I m (n = m QQ T = I n ) A Q q j r a j r q j Q T r = Q T (y Ax) = 0 (10.14) A QR Rx = Q T y (10.15) y q j y A m + 1 R

10. 4 (10.11) QR q 1 = q 2 = q 3 = r 33 = ε 1 [1, ε, 0, 0]T 1 + ε 2 1 (1 + ε2 )(2 + ε 2 ) [ε, 1, 1 + ε2, 0] T 1 (2 + ε2 )(3 + ε 2 ) [ε, 1, 1, 2 + ε2 ] T r 11 = 1 + ε 2 1 r 12 = r 13 = 1 + ε 2 2 + ε 2 r 22 = ε 1 + ε 2 ε r 23 = (1 + ε2 )(2 + ε 2 ) 3 + ε 2 2 + ε 2 ε r 33 r 22 ε x 3 x 2 ε ε 2 (10.1) (10.6) QR q 1 = 1 n [1 1 1] T q 2 = r 11 = n 1 (tj t) 2 [t 1 t t 2 t t n t ] T r 12 = n t r 22 = (tj t) 2 q T 1 y = n y q T [ 2 y q1 (q T 1 y) ] (tj t )(y j y) = (tj t ) 2 (10.4) A m + 1 y QR m + 1 q j a (m) m+1 r a (m) m+1 r 2 m+1,m+1 S (10.11) QR ( 6) n m A 0 (6.9) A [ ] P m P m 1 P 1 A = Q T R A = 0 R m m 0 (n m) m Q n n Q T Q = QQ T = I n Q y [ ] z = z (1) z (2) = Q T y z (1) m z (2) n m [ ] [ ] R z (1) x = Q T r (10.16) 0 z (2) Q r Q T r x n m (z (2) ) r m Rx = z (1) (10.17) r 2 = z (2) 2 Q R z A y

10. 5 r y y i σ 2 i ri 2 (10.7) σ 2 i S w = n r 2 i σ 2 i=1 i (10.18) C 1/2 y = σ 1 1 0 0 0 σ2 1........ 0 0 σ 1 n (10.19) (10.18) S w S w = (y Ax) T C 1 y (y Ax) (10.20) y = C 1/2 y y A = C 1/2 y A (10.21) S w S w = (y A x) T (y A x) (10.22) (10.7) A A y y y y y = E[y] y σy 2 = E[(y y) 2 ] y n y y C y = E[(y y)(y y) T ] (10.23) C y n n (C y ) ij = E[(y i y ) (y j y j )] y i y j C y y i y j (C y ) ij = 0 i j C y x (10.7) y (10.9) ˆx = By (10.24) B = (A T A) 1 A T ˆx E[ˆx] = E[By] = By ˆx C ˆx = E[(ˆx x)(ˆx x) T ] = E[B(y y)(y y) T B T ] = BE[(y y)(y y) T ]B T = BC y B T (10.25) y i σ0 2 C y = σ0i 2 n

10. 6 ˆx C ˆx = σ 2 0BB T = σ 2 0(A T A) 1 (10.26) (A T A) 1 A QR (A T A) 1 = (R T Q T QR) 1 = (R T R) 1 = R 1 R T R 1 R T R 1 R R R 1 R R ij R 1 R 1 ij j k=1 R ik R 1 kj = δ ij i j R 1 ij i = m, m 1,, 1 R 1 ii = 1 R ii R 1 ij = R 1 ii j k=i+1 j = i + 1, i + 2,, m R ik R 1 kj (10.27) i = m R 1 ij R ij (10.22) ˆx = B y B = (A T A ) 1 A T (10.28) = (A T C 1 y A) 1 A T Cy 1 A (10.21) C ˆx = B C y B T = (A T A ) 1 (10.29) A QR A = Q R C ˆx = R 1 R T (10.30) ˆx ŷ = Aˆx (10.31) Cŷ Cŷ = AC ˆx A T C ˆx (10.29) A = C 1/2 y A = Q R Cŷ = C 1/2 y Q Q T C 1/2 y (10.32) Q n m Q Q T C y σ 2 i Q q ij C ŷ (Cŷ) ii ŷ i (Cŷ) ii = σ 2 i m j=1 q 2 ij Q Q 1 1 (Q n n 1 ) (Cŷ) ij σ 2 i = (C y ) ii (10.33) ( ) ˆx R R ŷ Q Q

10. 7 C y (10.20) (10.21) C 1/2 y C y L C y = LL T (10.34) l 11 0 0. l 21 l.. 22. L =.... 0 l n1 l nn L ( 4) C y σ ij k = 1, 2,, n k 1 l kk = σ 2 kk j=1 l 2 kj l ik = 1 ( k 1 ) σik 2 l ij l kj l kk j=1 i = k + 1, k + 2,, n (10.35) C y C 1 y = (LL T ) 1 = L T L 1 (10.20) S w = L 1 (y Ax) 2 (10.21) y = L 1 y A = L 1 A C 1/2 y L 1 C 1/2 y L x µ σ 2 p(x µ, σ 2 ) = 1 [ exp (x ] µ)2 2πσ 2 2σ 2 (10.36) x 1 x 1 µ p(x 1 µ, σ 2 ) x 1 µ L(µ, σ 2 x 1 ) = p(x 1 µ, σ 2 ) (10.37) µ x 1 n x 1, x 2,, x n L(µ, σ 2 x 1, x 2,, x n ) = p(x 1 µ, σ 2 ) p(x 2 µ, σ 2 ) p(x n µ, σ 2 ) [ = (2πσ 2 ) n/2 exp 1 n 2σ 2 (x i µ) 2] (10.38) i=1 µ µ l(µ, σ x 1, x 2,, x n ) = log L(µ, σ 2 x 1, ) = n 2 log(2πσ2 ) 1 n 2σ 2 (x i µ) 2 (10.39) i=1 µ L n (x i µ) 2 = min i=1 µ 0 n (x i µ) = 0 i=1 ˆµ = 1 n n x i (10.40) i=1

10. 8 σ 2 (10.39) σ 2 0 n σ 2 + 1 σ 4 n i=1 (x i µ) 2 = 0 σ 2 ˆσ 2 = 1 n (x i ˆµ) 2 (10.41) n i=1 x i ˆµ ˆσ 2 E[ˆµ] = µ (10.42) (unbiased estimator) ˆσ 2 E[ˆσ 2 ] = n 1 n σ2 (10.43) σ 2 (10.41) ˆσ 2 σ 2 ˆσ 2 = 1 n 1 n (x i ˆµ) 2 (10.44) i=1 ˆx C ˆx (10.26) (10.29) x x = ˆx ( ˆx ) [ ] σ1 2 ρσ 1 σ 1 σ 2 C ˆx = (10.46) ρσ 1 σ 2 σ2 2 σ 2 1, σ 2 2 ˆx 1, ˆx 2 ρ ˆx 1 ˆx 2 C 1 ˆx = 1 1 ρ 2 δx i = x i ˆx i [ 1/σ 2 1 ρ/σ 1 σ 2 ρ/σ 1 σ 2 1/σ 2 δs = (x ˆx) T C 1 1 ˆx (x ˆx) = ( δx1 δx ) 2 2ρ + σ 1 σ 2 ( δx2 σ 2 ) 2 ] ] [( δx1 ) 2 1 ρ 2 σ 1 ˆx x 2 ˆµ E[(ˆµ µ) 2 ] = 1 n σ2 n 0 ˆµ = 0.95 x 1 = -0.90 x ˆx C ˆx l(x) = 1 2 (x ˆx)T C 1 ˆx (x ˆx) + (10.45) Fig. 10.1 Fig. 10.1 σ 2 /σ 1 = 2.5/2.0 δs = 1 ρ = 0.95( ) ρ = 0.9( )

10. 9 x 1 x 2 ρ = 0.95 δx 1 δx 2 QR R x 1, x 2, x 3 π x i e i 3 (x i + e i ) = π i=1 3 e i 2 = min i=1 e i e i x i + e i m x j n b i1 x 1 + b i2 x 2 + + b im x m = f i i = 1, 2,, n n < m Bx = f (10.47) B n m f n x x e x + e B(x + e) = f (10.48) n < m x e e 2 = min (10.49) (10.48) e 2 (10.48) Be = f Bx y (10.50) y x e T e + λ T (y Be) +(y Be) T λ = min (10.51) λ n e 0 δe T e + e T δe λ T Bδe δe T B T λ = 0 δe e = B T λ (10.51) λ (10.50) e (10.50) Be = BB T λ = y λ e = B T ( BB T ) 1 y (10.52) x ˆx = x + B T ( BB T ) 1 (f Bx) (10.53) x

10. 10 B = [1 1 1] f = π y = π x 1 x 2 x 3 BB T = 3 ˆx i = x i + 1 3 y (10.52) QR B T = QR (10.54) Q m n Q T Q = I n R n n (10.52) e = QR T y (10.55) (1982) :

11. 1 11 6 n m (n m) A A = UΛV T n m Ax = y S n n S T = S (11.1) S S T = S λ i u i Su i = λ i u i i = 1, 2,, n (11.2) u T j u i = 0 λ i λ j { u T 1 i = j j u i = δ ij = (11.3) 0 i j n n U = [u 1 u 2 u n ] (11.4) (11.3) U U T U = I n = UU T (11.5) I n n n U (11.2) SU = UΛ (11.6) Λ n n λ 1 0 λ 2 Λ =... 0 λ n (11.5) (11.7) S = UΛU T U T SU = Λ (11.8) S (11.8) 0 r 0 n r 0 S ( ) r λ i 0 λ i = 0 1 i r r + 1 i n (11.6) UΛ = [λ 1 u 1 λ 2 u 2 λ r u r 0 0] n r 0 U T U T n r U n r 0 UΛU T U n r U r = [u 1 u 2 u r ] S = U r Λ r U T r (11.9) Λ r 0 (11.7) n = r U r n r U T r U r = I r r < n U r U T r I n (r < n)

11. 2 S Sx = y (11.10) x y u i x = i y = i x iu i x i = u T i x y iu i y i = u T i y (11.11) (11.10) Sx = i n λ i x iu i = i=1 x isu i = i n y iu i i=1 x i λ i u i u i λ i x i = y i 1 i n (11.12) 0 x i = y i λ i x i (11.11) x n u T i x = y u i = UΛ 1 U T y (r = n) (11.13) λ i i=1 (11.8) (11.10) 0 (11.12) i > r 0 x i = y i r < i n y i = 0 (11.10) y i = 0 i > r (compatibility condition) u T i y = 0 i > r (11.14) S (11.10) y y u 1 u r 0 r u T i x = y u i + λ i i=1 n i=r+1 x iu i (11.15) x i S i > r u i S = [ 0 A A T 0 n m A ] (11.16) S (n + m) (n + m) Sw i = λ i w i (11.17) n + m w i n + m n u i m v i [ ] w i = u i v i (10.17) [ ] [ ] [ 0 A = λ i A T 0 u i v i u i Av i = λ i u i A T u i = λ i v i (11.18) u i v i A T Av i = λ 2 i v i AA T u i = λ 2 i u i (11.19) A T A m m m v i n ]

11. 3 (11.17) w i n + m (11.19) w i S (11.18) (λ i, u i, v i ) (λ i, u i, v i ) ( λ i, u i, v i ) (11.18) p (λ i, u i, v i ) λ i > 0 i = 1, 2,, p (11.20) ( λ i, u i, v i ) i = 1, 2,, p (11.21) 2p n + m 2p 0 0 u j v j (11.19) AA T u j = 0 j = p+1, p+2,, n A T Av j = 0 j = p+1, p+2,, m (0, u j, 0) j = p+1, p+2,, n (0, 0, v j ) j = p+1, p+2,, m (11.22) (11.20) (11.21) (11.22) p + p + (n p) + (m p) = n + m w i w T i w j = 0 i j w i w j (11.20) u T i u j + v T i v j = 0 i j w j (11.21) u T i u j v T i v j = 0 i j (11.20) (11.21) u i v i 0 w i u i v i u T i u j = δ ij v T i v j = δ ij (11.23) S S A u i v i U = [u 1 u 2 u p u p+1 u n ] V = [v 1 v 2 v p v p+1 v m ] (11.24) U n n V m m p 0 Av i = λ i u i AV = UΛ (11.25) Λ λ 1, λ 2,, λ p 0 n m U V V T U T A = UΛV T U T AV = Λ (11.26) 0 U p = [u 1 u 2 u p ] V p = [v 1 v 2 v p ] (11.27) p p λ 1 0 λ 2 Λ p =... 0 λ p (11.28) A = U p Λ p V T p (11.29) U p V p U T p U p = I p V T p V p = I p

11. 4 U p U T p = I n Λ 0 r = min(n, m) λ i 0 (1 i r) A p A ( ) A p 1 1 p r = min(n, m) p = r () p < r max λ i 1 i r κ(a) = min λ i 1 i r (11.30) A A T A AA T (11.19) κ(a T A) = κ(aa T ) = κ(a) 2 A 1 1 A = ε 0 0 ε A T A AA T A T A : λ 2 = 2 + ε 2, ε 2 AA T : λ 2 = 2 + ε 2, ε 2, 0 A p = 2 u 1 0 λ = 2 + ε 2, ε 2 0 ε U = ε 1 1 ε 1 1 [ ] 1 1 V = 1 1 n m A Ax = y (11.31) x y x m m v i n y u i x = m x iv i y = i=1 n y iu i (11.32) i=1 x i x v i y i y u i y y i y i = u T i y x (11.31) m m Ax = x iav i = x iλ i u i i=1 i=1 u i λ i x i = y i i = 1, 2,, p (11.33) 0 = y i i = p+1, p+2,, n (11.34) (11.31) (11.31) (11.33) x p (11.32) x x = p i=1 y i λ i v i = p i=1 u T i y v i λ i = V p Λ 1 p U T p y (11.35) (11.34) x i = y i /0 (i > p)

11. 5 y y u i (i > p) u T i y = 0 i = p+1, p+2,, n (11.36) x y (11.31) p n y Ax = (y i λ i x i)u i + y iu i i=1 i=p+1 x (11.31) (11.33) n y Ax = y iu i (11.37) i=p+1 (11.34) 0 (11.31) x i (i > p) A m x = x + x iv i (11.38) i=p+1 x i x (11.37) x x x n, m, p (n = m) n = m p = n = m ( ) 0 A (11.33) i = n (11.35) x = V Λ 1 U T y y (11.31) p < n (11.38) (11.37) (11.36) 0 x n > m n > m p = m (11.38) x (11.38) (11.36) 0 n < m 10 (10.50) p = n (11.37) 0 (11.38) 0 A (i) AXA = A (ii) XAX = X (11.39) (iii) (AX) T = AX (iv) (XA) T = XA X A A 1 4 (11.35) (11.31)

11. 6 (11.31) A = V p Λ 1 p U T p (11.40) A (11.39) 4 X = A A Λ p 0 0 q (11.40) A q = V q Λ 1 q U T q (11.41) A m n A A (11.29) A qa = (V q Λ 1 q U T q )(U p Λ p V T p ) q p U T q U p u i q q p q 0 U T q U p Λ p V T p = Λ q V T q (i) AA qa = U q Λ q V T q (ii) A qaa q = A q (11.42) (iii) (iv) AA q = U q U T q A qa = V q V T q q = p (11.29) (i) A A q < p (ii) (iv) (11.39) A q (11.31) x = A qy (11.43) x (11.31) x = R q x R q = A qa (11.44) x i x R q i x R q R q A q (11.42) R q = V q V T q q = m R q x y y = Ax (11.45) (11.43) y = S q y S q = AA q (11.46) S q y y x x S q y y A q (11.42) S q = U q U T q (11.47) q = n n > m p q S q (11.31) x Ax y y x (11.43)

11. 7 10 x C x C x = A qc y (A q) T (11.48) C y 10 y y i σ 2 i C y σ 2 i x i (C x ) ii = n ( q k=1 j=1 v ij u kj λ j ) 2σ 2 k (11.49) u ij, v ij u j, v j u j = [u 1j u 2j u nj ] T q 1 6 11 2 7 12 A = 3 8 13 4 9 14 5 10 15 λ = 35.127223 2.646397 5.15e-16 0 p = 2 q = p = 2 R 2 S 2 0.83 0.33 0.16 R 2 = 0.33 0.33 0.33 0.16 0.33 0.83 S 2 = 0.6 0.4 0.2 0.0 0.2 0.4 0.3 0.2 0.1 0.0 0.2 0.2 0.2 0.2 0.2 0.0 0.1 0.2 0.2 0.4 0.2 0.0 0.2 0.4 0.6 R 2 0.16 0.16666 0.42 0.36 0.315 0.28 10 1.26 0.84 0.63 0.504 0.84 0.63 0.504 0.42 A = 0.63 0.504 0.42 0.36 0.504 0.42 0.36 0.315 2.5582006 λ = 0.1799580 6.2085986e-3 9.9670848e-5 (p = 4) κ = 2.6 10 4 A T A 7 10 8 q = p = 4 S q S 4 = 0.99982 0.00211 0.00739 0.00986 0.00444 0.97465 0.08873 0.11831 0.05324 0.68944 0.41408 0.18634 0.44789 0.24845 0.88820 10 x = [1/2, 1/3, 1/4, 1/5] T

11. 8 y = Ax U V A 4 x 0.50000119 x = 0.33332490 0.25001625 0.19999081 Lanczos, C. (1961) : Linear Differential Operators, Van Nostrand.

12. 1 12 ( ) ( ) 1/100 30 Knuth TEX The Art of Computer Programming () TEX PC9801 MicroTEX (straight insertion) n a(1), a(2),, a(n) k a(1), a(2),, a(k) a(1) a(2) a(k) a(k+1) a(k) a(k+1) a(k+2) a(k) < a(k+1) a(k) a(k+1) a(k 1) a(k+1) a(k+1) do k=1, n-1 ak=a(k+1) do i=k, 1, -1 if( a(i).lt.ak ) a(i+1) = a(i) else goto 1 endif enddo 1 a(i+1)=ak enddo then k k/2 n 1 k=1 k 2 = 1 4 n(n 1) 1 4 n2 (12.1) n 100 1000 n 100 a(k+1) a(1) a(k) a(l) l = k + 1 2 l a(l) > a(k + 1) a(l) a(k) a(k+1) log 2 k n 1 log 2 k = log 2 (n 1)! k=1 log 2 (n 1)! n log 2 n n = 10 3. = 2 10 n 2 /4 25:1

12. 2 do k=1, n-1 ak=a(k+1) if( a(k).ge.ak ) goto 3 if( a(1).lt.ak ) then i=1 else il=1 ir=k 1 i=(il+ir+1)/2 if( i.eq.ir ) goto 2 if( a(i).gt.ak ) then il=i else ir=i endif goto 1 endif 2 do j=k, i, -1 a(j+1)=a(j) enddo a(i)=ak 3 enddo n 3/2 a(i) a(1) a(3) a(5) a(2) a(4) a(6) a(i) a(i + 1) a(1), a(2), a(3), n k a(1) a(n k +1) a(2n k +1) a(2) a(n k +2) a(2n k +2) a(n k ) a(2n k ) a(3n k ) n k 1 n k = 3k 1 2 n 1 = 1 n k+1 = 3n k + 1 (12.2) n > n k n k nk=1 do k=1, n nk=3*nk+1 if( nk.gt.n ) goto 1 enddo 1 nk=nk/3 if( nk.le.0 ) do i=nk+1, n ai=a(i) do return j=i, nk+1, -nk if( a(j-nk).lt.ai ) a(j)=a(j-nk) else goto 2 endif enddo 2 a(j)=ai enddo goto 1 then (heap) a(i/2) a(i) 1 i n (12.3) i/2 ( )

12. 3 1 2 3 4 5 6 7 n log 2 n (12.4) 8 9 10 11 12 13 14 15 a(1) a(2), a(3) () a(2) a(4), a(5) a(3) a(6), a(7) a(2) a(3) 2 2 a(1) a(n) a(1) (12.3) a(2) a(3) a(1) a(3) a(6) a(7) a(3) a(n) a(1) a(n 1) a(1) a(n) n = 15 a(14) a(15) a(7) a(7) a(12), a(13), a(6) a(6) (12.3) a(3) a(6), a(7) a(3) a(7) a(3) a(7) a(7) a(14), a(15) il=n/2+1 ir=n if( il.gt.1 ) then il=il-1 ai=a(il) else ai=a(ir) a(ir)=a(1) ir=ir-1 if( ir.eq.1 ) then a(1)=ai return endif endif i=il j=2*il 2 if( j.gt.ir ) goto 3 if( j.lt.ir.and. a(j).gt.a(j+1) ) * j=j+1 if( ai.gt.a(j) ) then a(i)=a(j) i=j j=2*j else j=ir+1 endif goto 2 3 a(i)=ai goto 1

12. 4 Fortran C C for do break Fortran Fortran x x a(i) x x a(j) x y a(i) y a(j) x (8, 10, 7, 9) (5, 3, 7) (1, 3, 0, 2) 10 x i 1 a(i) < x j n 1 a(j) > x a(i) a(j) i j x 2 log 2 n n = 10 6. = 2 20 40 Numerical Recipes in C ir=n il=1 is=0 1 if( ir-il.lt.m ) then do j=il+1, ir aj=a(j) do i=j-1, il, -1 if( a(i).ge.aj ) goto 2 a(i+1)=a(i) enddo i=il-1 2 a(i+1)=aj enddo if( is.eq.0 ) return ir= istck(is) il=istck(is-1) is=is-2 else k=(il+ir)/2 tmp=a(k) a(k)=a(il+1) a(il+1)=tmp if( a(il+1).lt.a(ir) ) then tmp=a(il+1) a(il+1)=a(ir) a(ir)=tmp endif if( a(il).lt.a(ir) ) then tmp=a(il) a(il)=a(ir) a(ir)=tmp endif if( a(il+1).lt.a(il) ) then tmp=a(il+1) a(il+1)=a(il) a(il)=tmp endif i=il+1 j=ir aj=a(il) 3 do k=il+1, ir if( a(i).lt.aj ) goto 4 i=i+1

12. 5 enddo 4 do k=il+1, ir if( a(j).gt.aj ) goto 5 j=j-1 enddo 5 if( j.ge.1 ) then tmp=a(i) a(i)=a(j) a(j)=tmp goto 3 endif a(il)=a(j) a(j)=aj is=is+2 if( is.gt.nstck ) then write(6,*) stack full return endif if( ir-i+1.ge.j-1 ) then istck(is)=ir istck(is-1)=i ir=j-1 else istck(is)=j-1 istck(is-1)=il il=i endif endif goto 1 istck nstck nstck 50 m m 10 m 1000 30 30 n m m n do k=1, n-1 ak=a(k+1) kx=min(k, m) do i=kx, 1, -1 if( a(i).lt.ak ) a(i+1)=a(i) else goto 1 endif enddo i=0 1 a(i+1)=ak enddo m then 1000 30 a(i) i 1 2 3 4 5 a(i) 45 68 93 70 50 indx(i) 3 4 2 5 1 indx(i) i 1 3 2 4 indx(i) = i a(i) a(i) indx(i) indx(i)

12. 6 do i=1, n indx(i)=i enddo nk=1 do k=1, n nk=3*nk+1 if( nk.gt.n ) goto 1 enddo 1 do k=1, n nk=nk/3 if( nk.le.0 ) return do i=nk+1, n ai=a(indx(i)) ia=indx(i) do j=i, nk+1, -nk if( a(indx(j-nk)).lt.ai ) * then indx(j)=indx(j-nk) else goto 2 endif enddo 2 indx(j)=ia enddo enddo n = 1000 n = 5000 n = 1000 n = 5000 Insertion 247,286 6,180,687 10 Shell 19,230 136,715 20 Heap 12,068 72,063 30 Quick 8,452 50,699 70 n 2 /4 n 3/2 n log 2 n n 1.1 n = 10, 000 PC 10 70

13. z 1 13 z α α x j x j k w k y j M y j = w k x j k (13.1) k= M w k x j y j w j = 1 (13.2) j= (13.1) y j w j x j x j (9.1) x j < x j j= 9 (13.1) x j = e ijω (9.1) ω 9 (13.1) y j = w k e i(j k)ω = e ijω w k e ikω k k = y j x j = k w k e ikω = W (ω) (13.3) W (ω) w j w j W (ω) W (ω) = W (ω) e iφ(ω) (13.4) y j = W (ω) e i(jω φ) W (ω) φ(ω) τ p (ω) = φ(ω) ω (13.5) (phase delay time) τ p t x j φ 0 w j φ 0 w j w j = w j (13.6) w j W (ω) = j= w j e ijω = w 0 + 2 w j cos jω j=1 φ 0 ( W (ω) < 0 ) 9 w j (13.6) (13.2) (13.2) W (0) = 1

13. z 2 Fig. 9.5 ω > 2π/M W M (ω) 0 M w j ( j M) h j = w j M j = 0, 1, 2,, 2M h j w j M h j H(ω) = = 2M j=0 h j e ijω = 2M j=0 w j M e ijω M w k e i(k+m)ω = W (ω)e imω k= M w j h j (13.5) τ p (ω) = M ω φ(ω) ω τ p (ω) ω y j y j = x j w k x j k k { 1 w j j = 0 v j = δ j w j = (13.7) w j j 0 y j = j v k x j k v j W (ω) V (ω) = 1 W (ω) x j y j (13.1) y j = k= h k x j k (13.8) δ j y j = k h k δ j k = h j h j h j < (13.9) j= h j H(ω) = j= h j e ijω (13.10) (13.10) h j = 1 2π π π H(ω)e ijω dω (13.11) (13.10) H(ω) (13.9) H(ω) h j 2 < j ( 9)

13. z 3 h j = 0 j < 0 (13.12) (causal) (physically realizable) x j ω ω/2 ω + ω/2 x j (13.11) x j (ω) = 1 2π ω+ ω/2 ω ω/2 X(σ) e i(jσ θ) dσ X(σ) x j θ(σ) x j X(σ). = X(ω) θ(σ). = θ(ω) + θ (ω)(σ ω) θ (ω) = dθ(ω)/dω x j (ω) =. ω 2π X(ω)e ijω sin(j θ ) ω/2 (j θ ) ω/2 sinc x j ω t g (ω) = dθ(ω) dω (13.13) (group arrival time) y j (13.8) H(ω) φ(ω) y j θ + φ x j φ(ω) (13.5) τ p (ω) y j (ω) j = θ + φ x j (ω) τ g (ω) = dφ(ω) (13.14) dω (group delay time) 0 (13.1) ω c { 1 ω < ω c H(ω) = (13.15) 0 ω > ω c (13.11) h j = 1 ωc e ijω dω = ω c sin jω c (13.16) 2π ω c π jω c 1/j ± (13.8) y j = M k= M h k x j k H(ω) H M (ω) = M j= M h j e ijω Fig. 13.1 t = 0.01s f c =10Hz (13.16) ω c = 2πf c t = 0.2π h j j = 0 j = 20 h j j 0 j = 20 (13.15) () h j w j w j h j w j (9.31) 6Hz

13. z 4 Impulse Response 5 10 15 20 (13.8) z y j = H(z)x j (13.19) H(z) z z = e iω (13.20) H(ω) H(ω) = H ( z = e iω) (13.21) 1.0 0 Frequency Response 10 20 Frequency z H z z H(z) z h j = 1 j dz H(z)z 2πi z z =1 (13.22) Fig. 13.1 z x j z z (13.18) H(z) zx j = x j 1 z A(z) B(z) z n x j = x j n (n ) (13.17) (13.9) H(z) = h k z k (13.18) k= z H(z) x j (13.17) H(z)x j = h k z k x j = h k x j k k k C(z) = A(z)B(z) (13.23) z A(z) a j B(z) b j C(z) c j C(z) = a j z j b k z k = a j b k z j+k j k j k = ( ) a j b n j z n n c n = j j a j b n j = k a n k b k (13.24) a j b j

13. z 5 z x n = cos nθe αn n 0, α > 0 (13.25) z X(z) = 1 2 = 1 ( 2 = n=0 ( e (iθ α)n + e (iθ+α)n) z n 1 1 e iθ α z + 1 1 e iθ α z 1 e α cos θz 1 2e α cos θz + e 2α z 2 (13.26) z z z A 1 (z) = z α A 2 (z) = 1 αz (13.27) α α A 1 (ω) 2 = A 2 (ω) 2 = 1 + α 2 2Re(αe iω ) A 1 (ω) = e iω α = A 1 (ω) e iφ1(ω) 1 da 1 A 1 (ω) dω = d dω ln A 1(ω) + i dφ 1(ω) dω dφ 1 /dω A 1 (z) 1 da 1 A 1 dω = i(1 αeiω ) A 1 (ω) 2 1 da 2 A 2 dω = i( α 2 αe iω ) A 2 (ω) 2 A 1 (ω) = A 2 (ω) dφ 1 (ω) dω dφ 2(ω) dω = 1 α 2 A 1 (ω) 2 ) α A 1 (z) α > 1 A 1 (z) dφ 1 (ω) dω < dφ 2(ω) dω α > 1 (13.28) A 1 (z) α A 1 (0) = A 2 (0) = 1 α φ 1 (0) = φ 2 (0) (13.28) φ 1 (ω) < φ 2 (ω) 0 < ω < π (13.29) A 1 (z) m z A(z) = a 0 + a 1 z + + a m z m (13.30) m α j A(z) A(z) = a m (z α 1 )(z α 2 ) (z α m ) A(z) z α j z α k (13.27) 1 α k z 2 m A(z) 0 z 1 (13.31) z m z z A(z) a 0 a m (wavelet) z (13.30) a j

13. z 6 α k α k A 1 (z) = (z α)(z α) = z 2 2Re α z + α 2 A 2 (z) = (1 αz)(1 αz) = 1 2Re α z + α 2 z 2 z α > 1 A 1 (z) A 2 (z) A 1 (0) = A 2 (0) φ 1 (ω) < φ 2 (ω) Fig. 13.2 7 (m = 7) 3 4 5 2 5 = 32 Fig. 13.2 16 16 Fig. 13.2 0 8 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 Fig. 13.2 7 (0 ) z B(z) B(z) z 0 ( ) B(z) n x j u j x j = B(z)u j u j x j u j u j = 1 B(z) x j (13.32)

13. z 7 1/B(z) (13.22) B(z) n n B(z) = b n (z β 1 )(z β 2 ) (z β n ) 1/B(z) 1 B(z) = n k=1 B k z β k (13.33) B k β j 1 1 = 1 1 z β j β j 1 β 1 j z = 1 ( ) 1 + β 1 j z + β 2 j z 2 + β j β j > 1 (z n n + ) β k < 1 1 = 1 1 z β k z 1 β k z 1 = z 1 + β k z 2 + β 2 kz 3 + β k < 1 (13.33) 1 B(z) = j= h j z j h j h j 0 h j (j < 0) j x j u j u j 1, u j 2, u j x j, x j 1, x j 2, k < 0 h k u j B(z) u j x k (k > j) B(z) B(z) B(z) ( ) B(z) 1/B(z) u j j < 0 0 x j j < 0 0 x j = B(z)u j n b 0 u j + b k u j k = x j j 0 k=1 x j u j B(z) u j = 1 n ) (x j b k u j k (13.34) b 0 k=1 j = 0, 1, 2, 8 B(z) x j u j B(z) u j u j x j u j B(z) u j u j u j u j (j )

13. z 8 (13.34) z A(z) = a 0 + a 1 z + a 2 z 2 + + a m z m B(z) = 1 + b 1 z + b 2 z 2 + + b n z n y j = A(z) B(z) x j (13.35) b 0 = 1 B(z) (13.34) m y j = a k x j k k=0 n b k y j k (13.36) k=1 j = 0, 1, 2, A(z)/B(z) x j y j = 1 1 z/2 x j h k = 2 k k 0 24 24 y j = x j + 0.5y j 1 (13.36) σ H(σ) 2 = 1 1 + (σ/σ c ) 2n (13.37) n n σ = 0 1 σ = ±σ c 1/2( ) σ ± 0 H(σ) n σ = ±σ c n H(σ) 2 { 1 σ < σ c 0 σ > σ c n (13.37) n = 1 (13.37) H(σ) 2 = H(σ) H(σ) = 1 [1 i(σ/σ c )][1 + i(σ/σ c )] 1 1 i(σ/σ c ) 1 1 + i(σ/σ c ) H(σ) h(t) = 1 2π H(σ)e iσt dσ h(t) H(σ) = 1 1 i(σ/σ c ) σ σ = iσ c R + R t < 0 e iσt 0 R + R + L -i c -plane Fig. 13.3 R -

13. z 9 H(σ) L+R + = 0 0 0 h(t) = 0 t < 0 t > 0 R = 2πi Res L+R Res H(σ) Res = iσ c t > 0 h(t) = σ c e σct t > 0 H(σ) = 1 1 + i(σ/σ c ) h(t) = 0 t > 0 n H(σ) 2 2n ( 2k + 1 ) σ/σ c = exp 2n iπ k = 0, 1, 2,, 2n 1 (13.38) z σ c n n H(σ) t < 0 h(t) = 0 H(σ) H(σ) H(σ) = n k=1 [ e iθ k i(σ/σ c )] 1 (13.39) θ k = 2k n 1 π 2n i (13.37) e iθ k (13.39) h(t) t H(σ) σ = 2 z 1 i t z + 1 (13.40) z σ t z z = 0 σ σ = 2i/ t z = σ = 2i/ t z z = e iω σ σ = 2 t tan ω 2 ω < π (13.41) ω σ. = ω t σ t ω

13. z 10 (13.40) (13.39) z z (13.39) z H(σ) σ π < ω < π σ = 0 z (13.40) (13.39) θ k θ n k+1 = θ k k n k + 1 z M (1 + z) 2 H(z) = k=1 b (k) 0 + b (k) 1 z + b(k) b (k) 0 = 1 + 4 b (k) 1 = 2 σ c t cos θ k + [ ( 2 ) 2 ] 1 σ c t b (k) 2 = 1 4 σ c t cos θ k + { n/2 n M = (n + 1)/2 n 2 z2 ( 2 ) 2 (13.42) σ c t ( 2 ) 2 σ c t n k = M 1 + z b (M) 0 + b (M) 1 z b (M) 0 = 1 + 2 σ c t b (M) 1 = 2 σ c t (13.43) (13.42) H(z) z n x j y j (13.36) (13.42) (13.42) H(z) = G M H k (z) (13.44) k=1 H k (z) = 1 + a(k) 1 z + a(k) 2 z2 1 + b (k) 1 z + b(k) 2 z2 1 y (1) j = H 1 (z)gx j y (2) j = H 2 (z)y (1) j y (M) j y (M) j y (k) j = y (k 1) j = H M (z)y (M 1) j y j + a (k) 1 y(k 1) j 1 + a (k) 2 y(k 1) j 2 b (k) 1 y(k) j 1 b(k) 2 y(k) j 2 (13.45) k = 1, 2,, M y (0) j = Gx j Fig. 13.4 6 Fig. 13.1 t = 0.01s f c = 10Hz f c 1/ 2 (13.42) σ c t (13.41) σ c t = 2 tan ω c 2 ω c = 2πf c t f = 10Hz H(ω) = 1/ 2 f = 8Hz 1 M = 3 13 40 Fig. 13.1 Fig. 13.4 t 0 6 f c 10 t

13. z 11 1.0 Amplitude 0 10 Butterworth Filter n = 6 10 20 Frequency H(z) y j m n z j = a k y j+k b k z j+k (13.46) k=0 j=1 j = N, N 1, N 2, j > N y j = 0 z Group Delay z j = H(z 1 )y j = H(z 1 )H(z)x j 0 10 20 Frequency Fig. 13.4 t (13.35) Fig. 13.4 10 H( ω)h(ω) = H(ω) 2 0 H(ω) H(ω) 2 (1978) 31, 240-263.

14. 1 14 ( 9) x(t) t t 0 1 E t [x(t)] lim T 2T + 1 T t= T x(t) = 0 (14.1) E t [ ] x(t) E t [ x(t) 2 ] 1 = lim T 2T + 1 T t= T x(t) 2 < (14.2) x(t) 0 x(t) r xx (τ) = E t [x(t + τ)x(t)] (14.3) x(t) τ r xx (0) (14.2) x(t) r xx ( τ) = r xx (τ) (14.4) x(t) r xx (τ) a k [ E t a k x(t + k) 2 ] 0 (14.5) k ] a k E t [x(t + k)x(t + l) k,l = a k r xx (k l)a l k,l a k a k r xx (k l)a l 0 (14.6) k,l r xx (τ) (positive definite) (non-negative definite) r xx (τ) (14.6) a k a 0 = 1 n a n = λ 0 (14.6) r xx (0) + λr xx (n) + λr xx ( n) + λ 2 r xx (0) 0 r xx ( n) = r xx (n) (14.4) λ n λ r xx (n) λr xx (n) = µ r xx (n) µ r xx (0) + 2µ r xx (n) + µ 2 r xx (0) 0 a l

14. 2 µ r xx (n) r xx (0) (14.7) 0 (14.5) x(t) t = 0 t = N x(t) r xx (τ) =. N τ 1 x(t + τ)x(t) τ 0 N τ + 1 t=0 t < 0 t > N x(t) = 0 r xx (τ) τ N + 1 r xx (τ) =. 1 N τ x(t + τ)x(t) τ 0 (14.8) N + 1 t=0 s(t) = e iω0t x(t) r ss (τ) = E t [e iω0(t+τ) e iω0t ] = e iω0τ E t [1] r ss (τ) = e iω0τ (14.9) ξ(t) 0 σ 2 ξ E t [ξ(t)] = 0 E t [ ξ(t) 2 ] = σ 2 ξ (14.10) ξ(t) E t [ξ(t + τ)ξ(t)] = 0 τ 0 (14.11) ξ(t) r ξξ (τ) = σ 2 ξδ τ (14.12) δ τ R xx (ω) = τ= r xx (τ) r xx (τ)e iωτ (14.13) x(t) r xx (τ) r xx (τ) = 1 2π π π τ = 0 r xx (0) = 1 2π π π R xx (ω)e iωτ dω (14.14) R(ω)dω r xx (0) R xx (ω) R xx (ω) x(t) R xx (ω) R xx (ω) ω (14.6) a k a k = e ikω 0 N k,l= N = N k N e ikω r xx (k l)e ilω N k+n k= N n=k N r xx (n)e inω

14. 3 (n = k l ) = 2N n= 2N r xx (n)e inω k Fig. 14.1 k { N n + N n 0 k = n N N n 0 k 1 = 2N + 1 n 2N n= 2N k ( 1 n ) r xx (n)e inω 0 2N + 1 N R xx (ω) () R xx (ω) 0 (14.15) 0 n 2N N N N k N 1 (14.9) δ(ω ω 0 ) 1 π δ(ω ω 0 )e iωτ dω = 1 2π π 2π e iω0τ (14.17) 1 2π e iω0τ e iωτ = δ(ω ω 0 ) (14.18) τ (14.17) R ss (ω) = 2πδ(ω ω 0 ) (14.19) x(t) y(t) r xy (τ) = E t [x(t + τ)y(t)] (14.20) x(t) y(t) x y E t [y(t + τ)x(t)] = E t [x(t τ)y(t)] r yx (τ) = r xy ( τ) (14.21) Fig. 14.1 k n (14.12) R ξξ (ω) = σξ 2 δ τ e iτω = σξ 2 (14.16) τ ω ( ) λ E t [ x(t + τ) + λy(t) 2 ] 0 (14.22) 0 r xx (0) + λr xy (τ) + λr xy ( τ) + λ 2 r yy (0)

14. 4 λr xy (τ) = µ r xy (τ) r xx (0) + 2µ r xy (τ) + µ 2 r yy (0) 0 (14.23) µ r xy (τ) 2 r xx (0)r yy (0) (14.24) (14.23) µ 0 (14.22) λ 0 (14.24) τ λ t x(t + τ) + λy(t) = 0 x(t) y(t) τ (14.7) (14.24) ρ xy (τ) = r xy (τ) rxx (0)r yy (0) (14.25) (14.24) ρ xy (τ) 1 (14.26) R xy (ω) = τ= r xy (τ) r xy (τ)e iωτ (14.27) x(t) y(t) E t [ N k= N [ x(t + k) + λy(t + k) ] e ikω 2 ] 0 N k,l= N [ r xx (k l) + λr xy (k l) ] +λr yx (k l) + λ 2 r yy (k l) e i(k l)ω 0 n = k l 2N n= 2N ( 1 n ) [rxx (n) 2N + 1 +2µr xy (n) + µ 2 r yy (n) ] e inω 0 N R xx (ω) + 2µR xy (ω) + µ 2 R yy (ω) 0 µ Rxy (ω) 2 R xx (ω)r yy (ω) (14.28) (14.7) (14.24) (14.28) γ xy (ω) = R xy (ω) Rxx (ω)r yy (ω) (14.29) γ xy (ω) 1 (14.30) x(t) h(t) y(t) y(t) = k h(k)x(t k) (14.31) r yx (τ) = E t [y(t + τ)x(t)] = k h(k)r xx (τ k)

14. 5 r oi (τ) = k h(k)r ii (τ k) (14.32) i o h(t) r ii (τ) R oi (ω) = H(ω)R ii (ω) (14.33) H(ω) [ r oo (τ) = E t h(k)x(t + τ k) ] h(l)x(t l) k l r oo (τ) = k,l = l h(k)r ii (τ + l k)h(l) r oi (τ + l)h(l) (14.34) r oi (τ) h(k) R oo (ω) = H(ω)R oi (ω) = H(ω) 2 R ii (ω) (14.35) (14.33) (14.35) γ oi (ω) = H(ω) H(ω) γ oi (ω) = 1 (14.36) 1 (14.28) x(t) y(t) AR (14.13) σξ 2 H(z) x(t) x(t) R xx (ω) (14.35) (14.16) R xx (ω) = σξ H(ω) 2 2 x(t) H(z) ξ(t) x(t) σ 2 ξ H(z) ξ H x Fig. 14.2 A(z) ξ(t) A(z) H(z) x A ξ Fig. 14.3 A(z) M a(k)x(t k) = x(t) + a(1)x(t 1) + k=0 +a(k)x(t k) + a(m)x(t M) = ξ(t) (14.37) a(0) = 1 a(0) = 1 a(k) ξ(t) x(t) = a(1)x(t 1) a(2)x(t 2) + ξ(t)

14. 6 x(t) x(t 1), x(t 2), ξ(t) y x y = α + βx + ε (14.37) (auto-regressive model) AR ξ(t) ξ(t) ξ(t) ξ(t) (14.34) M r ξξ (τ) = a(k)r xx (τ + l k)a(l) (14.38) k,l=0 0 M r ξξ (0) = a(k)r xx (l k)a(l) k,l=0 a(k) a(l) 0 M a(k)r xx (l k) = 0 l = 1, 2, 3,, M (14.39) k=0 M a(k) M r ξξ (0) l = 0 0 M σξ 2 = r ξξ (0) = a(k)r xx ( k) (14.40) k=0 (14.38) R ξξ (ω) = A(ω) 2 R xx (ω) A(ω) a(k) ξ(t) R ξξ (ω) = σ 2 ξ R xx (ω) = σ2 ξ A(ω) 2 (14.41) x(t) x(t) x(t) r xx (τ) (14.39) a(k) (14.40) σ 2 ξ a(k) σ ξ 2 (14.41) AR M (14.39) (14.39) (14.40) M a(k)r xx (l k) = σξδ 2 l k=0 l = 0, 1, 2,, M [a(0) a(1) a(m)]r M = [σ 2 ξ 0 0] (14.42) R M (M + 1) (M + 1) R M = r(0) r(1) r(m) r( 1) r(0) r(m 1)...... r( M) r( M +1) r(0) (14.43) r xx (k) r(k)

14. 7 r( k) = r(k) R M R T M = R M (14.44) R M (14.42) m [a m (0) a m (1) a m (m)]r m = [α m 0 0] a m (0) = 1 (14.45) a m (k) α m R m (14.43) (m + 1) (m + 1) m = 0 α 0 = r(0) (14.45) m + 1 m (14.45) 0 [a m (0) a m (1) a m (m) 0]R m+1 = [α m 0 0 γ m ] (14.46) m + 1 (14.45) m γ m = a m (k)r(m k + 1) (14.47) k=0 a m (k) γ m (14.46) [0 a m (m) a m (1) a m (0)]R m+1 = [γ m 0 0 α m ] (14.48) r( k) = r(k) (14.46) (14.45) 0 (14.48) u m (14.46) 0 γ m + u m α m = 0 (14.49) u m (14.45) m m + 1 a m+1 (k) = a m (k) + u m a m (m k + 1) k = 1 m a m+1 (m + 1) = u m (14.50) α m+1 = α m + u m γ m = α m (1 u m 2 ) m = 0 m = 0 m = 1, 2,, M 1 (14.42) (14.42) (14.45) α M σξ 2 α m α m ξ(t) (14.42) r(k) (14.42) α m α m R m 1 + δ r(0) 1 + δ δ 0.01 δ r(0) 1 (14.47) (14.49) (14.50) 2m + 2 m = 0 M 1 M(M + 1) (14.42) M 3 M AIC M M M

14. 8 ξ(t) M M N ξ(t) t 0 a(k) M AIC m AIC AIC m = N(ln 2πα m + 1) + 2(m + 1) (14.51) N x(t) α m m m AIC m m m AIC 1 2 AIC N α m N AIC x y z 0 1.0 2.0 (s) Fig. 4(a) x t t = 2msAR y t = A(z)x ty t AR z t = A(z)y t 3 Fig. 4(a) t = 2 ms Hz 20 30 Hz 2 2.4 1200 Fig. 4(b) r x 30 60 ms 100 ms 0.33(r(0) ) 140 ms M = 30 AR a(k) a(k) a(1) k 1 k > 10 M = 30 a(k) ξ(t) (14.41) R xx (ω) Fig. 4(c) 4 Hz 50 Hz (14.41) ξ(t) a(k) (14.37) ξ(t) Fig. 4(a) x(t) AR Fig. 14.4(b) r y 0

14. 9 10 ms a(k) a y a x ξ(t) Fig. 14.4(c) P y P x 3 50 P x 1.0e-5 1.0e-6 1.0e-7 1.0e-8 r x P y 1.0e-9 a x 20 40 60 (ms) 0 50 100 (Hz) r y a y 20 40 60 (ms) Fig. 14.4(c) P x M = 30 M = 6 P y ξ(t) m AIC m AIC 0 = N ln(α m /α 0 ) + 2m r z 20 40 60 (ms) Fig. 4(b) x(t) (r x) (a x)ξ(t) (r y) (a y)ξ(t) AR a x Hz ξ(t) Fig. 14.4(b) a y ξ(t) Fig. 14.4(a) z ξ(t) x Fig. 14.4(b) r z 3 0 M = 30 AIC x(t) α m AIC m AIC m α m AIC m α m m = 6 AIC m = 6 m = 9 m α m AIC m m α m AIC m 0 8.754e-8 0.0 10 2.574e-9-4211.8 1 3.388e-9-3900.3 11 2.573e-9-4210.3 2 3.306e-9-3927.6 12 2.573e-9-4208.4 3 2.883e-9-4089.9 13 2.573e-9-4206.4 4 2.656e-9-4186.1 14 2.573e-9-4204.5 5 2.593e-9-4213.0 15 2.572e-9-4203.0 6 2.587e-9-4213.7 16 2.569e-9-4202.2 7 2.587e-9-4212.0 17 2.567e-9-4201.5 8 2.582e-9-4212.2 18 2.566e-9-4199.9 9 2.577e-9-4212.5 19 2.566e-9-4197.9 AIC m m = 6 M = 6 x(t) Fig. 14.4(c) P x M = 30

14. 10 r x 6 M = 6 (14.37) z x(t) = 1 A(z) ξ(t) 1/A(z) a(k) (14.39) 13 1/A(z) h(k) ξ(t) x(t) h(k) r xx (τ) = σ 2 ξ h(τ + k)h(k) k=0 M r xx (τ) M = 30 1/A(z) Fig. 14.4(b) r x AR AR x(t) y(t) z(t) n x 1 (t), x 2 (t),, x n (t) x(t) = [x 1 (t) x 2 (t) x n (t)] T (14.52) r ij (τ) = E t [x i (t + τ)x j (t)] r(τ) = E t [x(t + τ)x (t)] r 11 (τ) r 12 (τ) r 1n (τ) r 21 (τ) r 22 (τ) r 2n (τ) =...... (14.53) r n1 (τ) r n2 (τ) r nn (τ) x = x T x x r(τ) x(t) (14.21) r(τ) r( τ) = r (τ) (14.54) (14.37) M ξ(t) = a(k)x(t k) a(0) = I n (14.55) k=0 x(t) n a(k) n n ξ(t) n I n n n r ξ (τ) = E t [ξ(t + τ)ξ (t)] [ ] = E t a(k)x(t + τ k)x (t l)a (l) k,l M = a(k)r(τ + l k)a (l) (14.56) k,l=0 r(τ) (14.53) (14.38) 1 0 ξ(t) σ 2 ξ = r ξ (0) M = a(k)r(l k)a (l) (14.57) k,l=0 1 a(k) σ 2 ξ σ 2 ξ a(k)

14. 11 a(k) σ 2 ξ (14.57) a(m) M δσ 2 ξ = δa(m) r(l m)a (l) + l=0 M a(k)r(m k)δa (m) k=0 (14.54) σ 2 ξ M a(k)r(l k) = 0 l = 1, 2,, M (14.58) k=0 n n a(0) = I n a(1) a(m) M a(k) (14.57) M σ 2 ξ = a(k)r( k) (14.59) k=0 (14.40) (14.56) P ξ (ω) r ξ (τ)e iωτ τ= = a(k)r(τ + l k)a (l)e iωτ k,l,τ = a(k)e iωk r(τ + l k)e iω(τ+l k) k τ [ a(l)e iωl l P ξ (ω) n n r ξ (τ) n n a(k) M A(ω) = a(k)e iωk (14.60) k=0 P ξ (ω) = A(ω)P (ω)a (ω) (14.61) P (ω) i, j r ij (τ) ξ(t) r ξ (τ) = 0 τ 0 P ξ (ω) = r ξ (0) = σ 2 ξ (14.61) P (ω) P (ω) = A 1 (ω)σ 2 ξ[a (ω)] 1 (14.62) r ij (τ) (14.58) a(k) (14.59) σ 2 ξ (14.60) (14.62) P (ω) P (ω) P (ω) (14.58) a(k) (14.59) σ 2 ξ [a(0) a(1) a(m)]r M = [σ 2 ξ 0 0] (14.63) [ ] M + 1 n n 0 n n R M (14.43) r(τ) n n r(τ) (14.44) (14.63) m [a m (0) a m (1) a m (m)]r m = [α m 0 0] [b m (m) b m (1) b m (0)]R m = [0 0 β m ] a m (0) = b m (0) = I n (14.64)

14. 12 a m (k), α m, b m (k), β m m = 0 α 0 = r(0) = β 0 a m (k) m α m b m (k) β m 1 b m (k) = a m (k) (14.64) [a m (0) a m (1) a m (m) 0]R m+1 = [α m 0 0 γ m ] (14.65) [0 b m (m) b m (1) b m (0)]R m+1 = [δ m 0 0 β m ] (14.66) R m+1 (14.65) m + 1 (14.66) m + 1 (14.64) m γ m = a m (k)r(m k + 1) k=0 m δ m = b m (k)r(k m 1) (14.67) k=0 (14.65) (14.66) (14.64) a m (k) (14.66) u m (14.65) γ m + u m β 0 (14.64) γ m + u m β m = 0 (14.68) u m u m n n u m u m (14.64) m m + 1 m a m+1 (k) a m+1 (k) = a m (k) + u m b m (m k + 1) k = 1 m a m+1 (m + 1) = u m (14.69) α m+1 = α m + u m δ m a m+1 (k) b m+1 (k) (14.65) n n v m (14.66) δ m + v m α m = 0 (14.70) v m b m+1 (k) b m+1 (k) = b m (k) + v m a m (m + 1 k) k = 1 m b m+1 (m + 1) = v m (14.71) β m+1 = β m + v m γ m m = 0 m = 1, 2,, M 1 α m β m α m = α m β m = β m (14.72) α m β m γ m δ m (14.69) (14.71) u m v m α m+1 = α m γ m β 1 m δ m β m+1 = β m δ m α 1 m γ m α m+1 α m+1 = α m δ m(β 1 m ) γ m α m β m α m+1 δ m = γ m (14.73)

14. 13 (14.72) (14.73) (AIC) AIC m = N(n log 2π + log α m + n) +n(n + 1) + 2n 2 m (14.74) N x(t) n α m α m 1 α m

15. 1 15 1 dy = f(x, y) (15.1) dx y(x n ) = y n x n+1 = x n + h y n+1 = y(x n+1 ) f n = f(x n, y n ) dy 1 dx = f 1(x, y 1, y 2 ) dy 2 dx = f 2(x, y 1, y 2 ) (15.2) [ ] [ ] y = y 1 y 2 f = (15.2) f 1 f 2 (15.3) dy = f(x, y) (15.4) dx (15.1) d 2 y dx 2 = f(x, y, y ) y 1 = y y 2 = dy dx [ d dx y 1 y 2 ] = f(x, y) = [ y 2 f(x, y 1, y 2 ) n n 1 ] (15.4) (15.1) y(x n + h) ( ) dy y(x n + h) = y n + h + 1 ( d 2 ) y dx n 2 h2 dx 2 n + 1 ( d 3 ) y 3! h3 dx 3 + (15.5) n n x n (15.1) dy = f(x, y) dx d 2 y dx 2 = df dx = f dy x + f y dx = f x + f y f (15.6) d 3 y dx 3 = f xx + 2f xy f + f x f y + fy 2 f + f yy f 2 f x, y (15.5) x n+1 y y n+1 f(x, y) (15.5) y n+1 = y n + hf(x n, y n ) (15.7) - y n = y(x n ) x n+1 y(x) y n+1 y n+1 = y n + γ 1 k 1 + γ 2 k 2 (15.8) k 1 = hf(x n, y n ) k 2 = hf(x n + αh, y n + βk 1 ) y n+1 (15.5) (15.6) y n+1 = y n + hf + 1 2 h2 (f x + f y f) + O(h 3 )

15. 2 f (x n, y n ) k 2 k 2 = hf + h 2 (αf x + βf y f) + O(h 3 ) (15.8) h 2 γ 1 + γ 2 = 1 γ 2 α = 1 2 γ 2 β = 1 2 α, β, γ 1, γ 2 (15.9) (15.6) h 3 (15.8) 2 (15.9) f(x, y) 2 h 2 2 2 (15.9) α = β β 3 α ( y n+1 = y n + 1 1 ) k 1 + 1 2α 2α k 2 (15.10) k 1 = hf(x n, y n ) k 2 = hf(x n + αh, y n + αk 1 ) α = 1 (Heun) α = 1/2 2 p y n+1 = y n + p γ i k i i=1 k i = hf ( i 1 ) x n + α i h, y n + β ij k j i = 1, 2,, p j=1 (15.11) α i, β ij, γ i (15.6) p - 4 (p = 4) 4 2 - y n+1 = y n + 1 6 (k 1 + 2k 2 + 2k 3 + k 4 ) + O(h 5 ) k 1 = hf(x n, y n ) k 2 = hf ( x n + 1 2 h, y n + 1 2 k ) 1 k 3 = hf ( x n + 1 2 h, y n + 1 2 k ) 2 k 4 = hf ( ) x n + h, y n + k 3 (15.12) f(x, y) y k x x n, x n + h x n + h/2 x f(x, y) h/2 h/2 h 1 /2, h 1 /2, h 2 /2, h 2 /2, f(x, y) y (15.1) dy dx = f(x) x y(x) = y n + f(x)dx x n y n+1 = y n + h 6 [ f(xn ) + 4f(x n +h/2) + f(x n+1 ) ] (15.12) k 2 = k 3 - (15.4) y f k i

15. 3 (15.12) k 4 h 4 - - - k 1 = hf(x n, y n ) z 1 = y n + ( 1 2 k 1 q 0 ) q 1 = q 0 + 3 ( 1 2 k ) 1 1 q 0 2 k 1 k 2 = h ( x n + 1 2 h, z ) 1 z 2 = z 1 + γ 2 (k 2 q 1 ) q 2 = q 1 + 3γ 2 (k 2 q 1 ) γ 2 k 2 k 3 = hf ( x n + 1 2 h, z ) 2 z 3 = z 2 + γ 3 (k 3 q 2 ) q 3 = q 2 + 3γ 3 (k 3 q 2 ) γ 3 k 3 k 4 = hf(x n + h, z 3 ) y n+1 = z 3 + 1 (1 3 2 k ) 4 q 3 q 4 = q 3 + ( 1 2 k ) 1 4 q 3 2 k 4 γ 2 = 1 1 1 = 2 2 + 2 = 0.292893 γ 3 = 1 + 1 = 2 + 2 = 1.70710 2 2 (15.13) q i q 0 0 q 4 q 0 q 1 = k 1 2q 0 -- - 4 5 h 1 y 1 (x + h) = y(x + h) + ch 5 + O(h 6 ) y(x + h) c x x + 2h 2h 1 y 1 (x + 2h) h 2 y 2 (x + 2h) y 1 (x + 2h) = y(x + 2h) + c(2h) 5 + O(h 6 ) y 2 (x + 2h) = y(x + 2h) + 2ch 5 + O(h 6 ) 2 2 h 5 y(x + 2h) = y 2 (x + 2h) + δ 15 + O(h6 ) (15.14) δ = y 2 (x + 2h) y 1 (x + 2h) 4 δ - (15.12) - - (15.13) 4 (p = 4)4 (15.11) p = 6 6 y n+1 = y n + 6 γ i k i + O(h 6 ) i=1 k i = hf ( i 1 ) x n + α i h, y n + β ij k j i = 1, 2,, 6 j=1 (15.15) 5 p 4 p p 1 p 2 6 5 5 4 4 4 y n+1 = y n + 6 γi k i + O(h 5 ) (15.16) i=1

15. 4 γi 0 4 yn+1 δ = y n+1 yn+1 = (γ i γi )k i (15.17) i γ i γ i γ i α i β ij γ i γ i γi 0 16/135 1/360 1/4 1/4 0 0 3/8 8/32 9/32 6656/12825 128/4275 12/13 1932/2197 7200/2197 7296/2197 28561/56430 2197/75240 1 439/216 8 3680/513 845/4104 9/50 1/50 1/2 8/27 2 3544/2365 1859/4104 11/40 2/55 2/55 Table 15.1 Fehlberg f(x, y) x f(x, y) Fig. 15.1 6 Cash-Karp 0 1 Fig. 15.1 α i Fehlberg Cash-Karp α i β ij γ i γ i γi 0 37/378 277/64512 1/5 1/5 0 0 3/10 3/40 9/40 250/621 6925/370944 3/5 3/10 9/10 6/5 125/594 6925/202752 1 11/54 5/2 70/27 35/27 0 277/14336 7/8 1631/55296 175/512 575/13824 44275/110592 253/4096 512/1771 277/7084 Table 15.2 Cash-Karp (15.14) (15.17) δ h 5 y n+1 δ 0 δ = ch 5 (15.18) δ δ 0 y n+1 y n+1 x n+1 x n+2

15. 5 h h 0 δ 0 δ 0 = ch 5 0 (15.18) c h 0 ( δ0 ) 1/5 h 0 = h δ δ 0 (15.19) δ δ > δ 0 y n+1 x n h 1 y(x n + h 1 ) δ 1 = ch 5 1 δ 0 h 1 δ 0 x n+1 h 1 x n x n + h h/h 1 x n + h δ 1 (h/h 1 ) δ 0 h 1 h h 1 δ 1 = chh 5 1 = δ 0 c h 1 = h ( δ0 δ ) 1/4 δ > δ 0 (15.20) (15.19) h 0 (15.20) h 1 y(x n + h) 0.9 - y(x n + h) y(x n ) (self-starting) h y n y n 1, y n 2, f(x, y) f n 1, f n 2, (15.1) y n+1 = y n + xn+1 x n f(x, y)dx f(x, y) h f n, f n 1, - ( ) y n+1 4 - y n+1 = y n + h 24( 55fn 59f n 1 +37f n 2 9f n 3 ) + O(h 5 ) (15.21) y n+1 - y n+1 = y n + h 24( 9fn+1 + 19f n 5f n 1 + f n 2 ) + O(h 5 ) (15.22) f n+1 f(x n+1, y n+1 ) y n+1 (implicit) - P : y n+1 = y n + h 24( 55fn 59f n 1 +37f n 2 9f n 3 ) (15.23) E : f n+1 = f(x n+1, y n+1) (15.24) C : y n+1 = y n + h 24( 9f n+1 + 19f n 5f n 1 + f n 2 ) (15.25) E : f n+1 = f(x n+1, y n+1) (15.26) P (predictor) C (corrector) E (evaluation) - y n+1 y n+1 f n+1 f n+1

15. 6 - y n+1 f n+1 (15.24) (15.25) P(EC)(EC) E - (15.23) (15.26) PECE y(x) x x + H z 0 = y(x) h = H n z 1 = z 0 + hf(x, z 0 ) z m+1 = z m 1 + 2hf(x + mh, z m ) (15.27) m = 1, 2,, n 1 y(x + H). = y n = 1 2 [z n + z n 1 + hf(x + H, z n )] h h y(x + H) H h (15.1) dy(x) dx = y(x + h) y(x h) 2h + O(h 2 ) 2 H y n y(x + H) = c i h 2i (15.28) i=1 h n y n/2 y n 4y n y n/2 3 h 4 n 2 H h f(x, y) y h n n = 2, 4, 6, 8, n h 2 h 2 = 0 5 Implicit Real*8 (a-h,o-z) Dimension x(50),y(50) x0 = 0 y0 = 1 yp0 = f(x0,y0) Do nn=1,10 n = nn*2 h = (xn-x0)/n x(nn) = h*h h2 = h*2 z1 = y0 z2 = z1+h*yp0 Do m=1, n-1 z0 = z1 z1 = z2 z2 = z0+h2*f(x0+m*h,z1) Enddo y(nn) = (z2+z1+h*f(xn,z2))/2 If( nn.gt.1 ) then Do i=nn-1,1,-1 y(i) = y(i+1)+(y(i+1)-y(i)) * *(x(nn)/(x(i)-x(nn))) Enddo Endif Write(6,6) n,(y(i),i=1,nn) 6 Format(i5,1p5e18.10/(5x,1p5e18.10)) Enddo Stop End Function f(x,y)

15. 7 Real*8 f,x,y f=-x*y Return End dy = xy y(0) = 1 dx H = 1 1 Do i=nn-1,1,-1 y(1) h 2 = 0 y(x) = e x2 /2 y(1) = e 1/2 = 0.606530659712 n y(1) y n RKG 2 0.6250 0.6250 0.606494 4 0.6035156 0.60888 0.606531 6 0.6065243 0.60738 0.606531 8 0.6065285 0.606975 0.6065308 10 0.606530639 0.606803 0.60653072 12 0.60653065892 0.606716 0.60653069 14 0.60653065970 0.6066649 0.606530678 16 0.60653065971 0.6066325 0.606530670 y(1) y n (15.27) h - - x = 1 (RKG) n = 14 10 y n y(1) 2 1 d 2 y = f(x, y) (15.29) dx2 y(x 0 ) = y 0 y (x 0 ) = y 0 y (x) = dy(x)/dy 2 d 2 y(x) y(x + h) 2y(x) + y(x h) dx 2 = h 2 + O(h 2 ) z 0 = y 0 z 1 = y 0 + h[y 0 + 1 2 hf(x 0, y 0 )] z m+1 2z m + z m 1 = h 2 f(x 0 + mh, z m ) y n = z n m = 1, 2,, n 1 (15.30) y n = 1 h (z m z m 1 ) + 1 2 hf(x 0 + H, z n ) y n y(x 0 + H) y n x 0 + H 1 2 m = z m+1 z m 1 0 = h[y 0 + 1 2 hf(x 0, y 0 )] z 1 = y 0 + 0 m = m 1 + h 2 f(x 0 + mh, z m ) z m+1 = z m + m (15.31) m = 1, 2,, n 1 y n (15.28)

16. 1 16 3 3 y = x 2 + c ( ) x 2 y 2 = c ( ) x 2 + y 2 = c ( ) 1 u t = D 2 u x 2 (16.1) D t 0 u(x, t 0 ) u(x, t) (16.1) 2 u t 2 = c2 2 u x 2 (16.2) c t 0 u(x, t 0 ) u/ t (16.2) 2 u x 2 + 2 u = ρ(x, y) (16.3) y2 (x, y) u(x, y) x x t u(x, t) u n j = u(j x, n t) 0 j J u n 0, u n J (16.1) ( ) n u t j = un+1 j u n j + O( t) (16.4) t ( 2 ) n u x 2 = un j+1 2un j + un j 1 ( x) 2 + O( x 2 ) (16.5) j (16.1) u n+1 j = u n j + β(u n j+1 2u n j + u n j 1) (16.6) 0 < j < J β = D t ( x) 2 (16.7) n n + 1 t (16.6) t u n j = λ n e ikj x (16.8) k x (16.6) λ = 1 4β sin 2 k x 2 (16.9) λ 1 1 λ 1 (16.9) k 0 β = D t ( x) 2 1 2 (16.10) t x (16.5) x (16.10) t 4 1 t = 0 t > 0 (16.1)

16. 2 (16.10) l τ = l2 D τ t = l2 D t = ( x)2 2D t 2l 2 ( x) 2 (16.10) τ ( l t 2 x ) 2 l x l l x τ/ t t (16.1) (16.6) (16.1) n (16.4) n n + 1 n + 1/2 ( ) n+1/2 u t j = un+1 j u n j + O( t 2 ) (16.11) t (16.4) n n + 1 2 u. 1 [ = θ(u n+1 x 2 ( x) 2 j+1 2un+1 j + u n+1 j 1 ) +(1 θ)(u n j+1 2u n j + u n j 1) ] (16.12) θ 0 θ 1 θ = 0 (16.6) u n+1 j = u n j + β [ θ ( u n+1 j+1 2un+1 j + u n+1 ) j 1 +(1 θ) ( u n j+1 2u n j + u n )] j 1 (16.13) u n+1 j n + 1 u n+1 j βθu n+1 j+1 + (1 + 2βθ)un+1 j βθu n+1 j 1 = β(1 θ) ( u n j+1 2u n j + u n ) j 1 (16.14) n+1 j+1, j, j 1 u n+1 j, j = 1, 2,, J 1 4 (16.8) (16.13) λ = 1 4β[λθ + (1 θ)] sin 2 k x 2 λ = 1 4β(1 θ) sin2 k x/2 1 + 4βθ sin 2 k x/2 λ 1 (16.15) 0 β sin 2 k x 2 2β(1 2θ) sin 2 k x 1 (16.16) 2 β 0 (16.10) θ θ < 1/2 β 1 2(1 2θ) 0 θ < 1 2 (16.17) β (16.10) θ θ = 1/2 (16.16) θ 1/2 θ 1/2 (16.14) θ = 1/2 - t

16. 3 - v t + v v x = D 2 v x 2 (16.18) u u = v (16.19) t x v (16.4) n u n j = v un j+1 un j 1 (16.20) t 2 x (16.8) u n+1 j λ = 1 i v t sin k x x t = 0 λ > 1 (16.20) (upwind differencing) v (16.20) u n+1 j u n v un j u j 1 v > 0 j = x t v un j+1 (16.21) un j v < 0 x λ = 1 2v t x i 2v t x λ 2 = 1 4 v t x sin k x 2 sin2 k x 2 ( 1 cos k x 2 v t x ) sin 2 k x 2 (16.21) v t x < 1 (16.22) (16.18) (ADI ) ( u 2 ) t = D u x 2 + 2 u y 2 (16.23) u(x, y, t) x y x u n j,l = u(j x, l x, n t) x y (16.5) (16.10) 0 D t ( x) 2 1 4 (16.24) - δ 2 xu n j,l = u n j+1,l 2u n j,l + u n j 1,l δ 2 yu n j,l = u n j,l+1 2u n j,l + u n j,l 1 - u n+1 j,l = u n j,l + 1 2 β[ δxu 2 n+1 j,l + δxu 2 n j,l +δyu 2 n+1 j,l + δyu 2 n ] j,l (16.25) x y (16.25) ADI (Alternating- Direction Implicit method,

16. 4 ) t u n+1/2 j,l u n+1 j,l = u n j,l + 1 2 β( δ 2 xu n+1/2 j,l = u n+1/2 j,l + δyu 2 n ) j,l + 1 2 β( δxu 2 n+1/2 j,l + δyu 2 n+1 ) j,l (16.26) u n+1/2 j,l l u n+1/2 j+1,l, un+1/2 j,l, u n+1/2 j 1,l u n+1 j,l+1, un+1 j,l, u n+1 j,l 1 - (16.18) - { u n+1/2 j = u n j v t (u n j un j 1 ) v > 0 x (u n j+1 un j ) v < 0 u n+1 j = u n+1 j + 1 2 β( δ 2 u n+1 j + δ 2 u n+1/2 j δ 2 u n j = u n j+1 2u n j + u n j 1 (16.27) n n+1/2 t t (16.22) 2 t t t ADI (16.26) t/2 t/2 x y t/2 (16.2) (16.5) u n+1 j 2u n j + u n 1 j = α 2( u n j+1 2u n j + u n ) j 1 α = c t x ) (16.28) (16.29) n u n+1 j u n j u n j = e i(kj x ωn t) (16.30) ω k (16.28) sin 2 ω t 2 sin ω t 2 = α 2 sin 2 k x 2 = α sin k x 2 (16.31) (16.31) k ω α 1 (16.31) sin ω t 2 > 1 (16.31) e iω t > 1 (16.28) α = c t x 1 (16.32) (16.28) 1 1 1 c t c t > x (16.2) c (16.28) ĉ (16.31) ĉ = ω k = 2 k t sin 1( α sin k x ) 2

16. 5 ĉ ( k x ) ( 1 c = 2 α sin 1 α sin k x ) 2 (16.33) ĉ k (16.2) (16.28) Fig. 16.1 (16.33) λ x x/λ ĉ/c α (j, l) u j,l u j,l 5 0 J = L = 4 1 9 l 7 8 9 1.0 0.95 0.9 0.8 4 5 6 1 2 3 j Relative Dispersion Fig. 0.9 0.8 0.7 0.6 0.7 0.5 0.1 0.5 0 0.1 0.2 0.3 0.4 0.5 x/ 16.1 λ ĉ/c α = c t/ x (16.3) 0 x J x 0 y L x u(x, y) (16.5) (16.3) u j+1,l + u j 1,l + u j,l+1 + u j,l 1 4u j,l = ( x) 2 ρ j,l (16.34) 0 < j < J 0 < l < L Fig. 16.2 u j,l x (16.34) Ax = b (16.35) b ρ j,l A A = 4 1 0 1 1 4 1 0 1 0 1 4 0 0 1 1 0 0 4 1 0 1 1 0 1 4 1 0 1 1 0 1 4 0 0 1 1 0 0 4 1 0 1 0 1 4 1 1 0 1 4 (16.36) 0 (16.34) (16.34) 4 0 1

16. 6 1 A A = E + F (16.37) E Ex = c F E E A (16.35) Ex = F x + b x (0) Ex (k) = F x (k 1) + b k = 1, 2, (16.38) x (k) k e (k) = x (k) x (16.39) (16.38) Ee (k) = F e (k 1) B = E 1 F e (k) = Be (k 1) (16.40) e (k) = B k e (0) e (k) k B λ i λ i < 1 ( 6) B λ i ρ(b) = max λ i (16.41) i B (16.38) e (k) ρ(b) (16.3) u t = 2 u x 2 + 2 u y 2 ρ (16.42) u(x, y, t) = 0 ( x, y ) t u n+1 j,l = u n j,l + t ( x) 2 [ u n j+1,l + u n j 1,l +u n j,l+1 + u n j,l 1 4u n j,l ( x) 2 ρ j,l ] t/( x) 2 (16.24) t t/( x) 2 = 1/4 u n+1 j,l = 1 [ u n j+1,l + u n j 1,l + u n j,l+1 4 ] +u n j,l 1 ( x) 2 ρ j,l (16.43) (relaxation method) (16.42) t ADI A A = D + L + U (16.44) D U L E = D (16.38) Dx (k) = (L + U)x (k 1) + b (16.45)

16. 7 D (16.34) u (k) j,l = 1 [ u (k 1) j+1,l 4 + u(k 1) j 1,l + u(k 1) j,l+1 ] +u (k 1) j,l 1 ( x)2 ρ j,l (16.46) 4 u (k) j,l (16.43) (j, l) u (k) j,l u (k) j,l u (k) u (k) j,l u (k 1) j,l Fig. 16.2 u (k) j,l = 1 [ u (k 1) j+1,l 4 + u(k) j 1,l + u(k 1) j,l+1 + u(k) j,l 1 ] ( x) 2 ρ j,l (16.47) - (D + L)x (k) = Ux (k 1) + b (16.48) (16.45) B = D 1 (L + U) ρ J ρ J = 1 ( cos π 2 J + cos π ) (16.49) L - ρ GS ρ GS = ρ 2 J J = L = 4 9 9 A (16.49) ρ J = 1 2 2 J (16.49) ρ J = cos π J 1 π2 2J 2 ρ J 1 1/e = 0.367 N e ρ Ne J = e 1 N e = 1 ln ρ J J ln ρ J ln ( 1 π2 2J 2 ) π 2 2J 2 N e 2 π 2 J 2 O(J 2 ) - J SOR - (16.48) (D + L)x (k) = (D + L)x (k 1) (D + L + U)x (k 1) + b = (D + L)x (k 1) Ax (k 1) + b x (k) r (k) = Ax (k) b (16.50) x (k) = x (k 1) (D + L) 1 r (k 1) (16.51) ω x (k) = x (k 1) ω(d + L) 1 r (k 1) (16.52)

16. 8 ω overrelaxation parameter SOR(successive overrelaxation) 0 < ω < 2 (16.34) 1 < ω < 2 (16.51) (16.34) r j,l = u j+1,l + u j 1,l + u j,l+1 + u j,l 1 4u j,l ( x) 2 ρ j,l (16.53) (16.51) u new j,l = u j,l + ω 4 r j,l (16.54) ω = 1 (16.47) SOR ω ρ J ω 2 ω opt = 1 + 1 ρ 2 J (16.55) ω SOR ρ SOR = ( ρ J 1 + 1 ρ 2 J ) 2 (16.56) J J ω opt 2 1 + π/j ρ SOR 1 2π J SOR N e N e 1 2π J SOR - J

17. 1 17 (Conjugate Gradient method, CG ) CG n n S S S T = S n a a T Sa > 0 (17.1) S Sx = c (17.2) F (x) = 1 2 xt Sx 1 2 (xt c + c T x) (17.3) x x T c F (x) (17.1) x F (x) F (x) (L 2 ) F (x) δf (x) = F (x + δx) F (x) = 1 2 δxt Sδx + 1 (δx T (Sx c) + [ δx T (Sx c) ] ) T 2 S F (x) Sx = c (17.2) F (x) x 0 x 1, x 2, x k p k x k+1 x k+1 = x k + α k p k (17.4) α k F (x k+1 ) F (x k + α k p k ) α k 0 F (x k + α k p k ) = α k p T k Sp α k k + 1 ( p T k (Sx k c) + [ p T k (Sx k c) ] ) T = 0 2 s k = c Sx k (17.5) (17.2) x k α k pt k s k α k = p T k Sp k (17.6) α k s k+1 (17.5) s k+1 = s k α k Sp k (17.7) (17.6) α k p T k s k+1 = 0 (17.8) x k+1 p k F (x k+1 ) F (x k )

17. 2 F (x) x k F (x) F (x k ) = c Sx k = s k (17.9) s k p k F (x) p k s k p k 1 p k = s k + β k 1 p k 1 (17.10) p k p k 1 S (conjugate) p T k Sp k 1 = p T k 1Sp k = 0 (17.11) (17.10) p T k 1S(s k + β k 1 p k 1 ) = 0 β k = pt k Ss k+1 p T k Sp k (17.12) CG CGS CGS (a) s 0 = c Sx 0 p 0 = s 0 (b) α k = s k 2 p T k Sp k (c) x k+1 = x k + α k p k (d) s k+1 = s k α k Sp k (17.13) s k+1 c (e) β k = s k+1 2 s k 2 (f) p k+1 = s k+1 + β k p k (b) α k β k (17.11) β k p T j s k = 0 j < k (17.14a) s T j s k = 0 j k (17.14b) p T j Sp k = 0 j k (17.14c) k j < k k = 1 p 1, s 1 s k+1 (17.13)(c) x k+1 p i k x k+1 = x j+1 + α i p i i=j+1 s k+1 k s k+1 = c Sx k+1 = s j+1 α i Sp i i=j+1 k p T j s k+1 = p T j s j+1 α i p T j Sp i i=j+1 (17.8) 0 j < i k 0 p T j s k+1 = 0 j < k j = k (17.8) (17.14a) (17.13)(f) 0 = p T j s k+1 = (s j + β j 1 p j 1 ) T s k+1 0 s T j s k+1 = 0 j k (17.4b) (17.13)(f) p T j Sp k+1 = p T j S(s k+1 + β k p k ) 0 p j (17.13)(c) p T j Sp k+1 = 1 ( ) T xj+1 x j Ssk+1 α j

17. 3 s j = c Sx j ( xj+1 x j ) T S = ( sj+1 s j ) T p T j Sp k+1 = 1 α j ( sj+1 s j ) T sk+1 = 0 s j p T j Sp k+1 = 0 j k (17.14c) α j 0 α j 0 (17.6) p T j s j = 0 (17.13)(f) 0 = p T j s j = ( s j + β j 1 p j 1 ) T sj = s j 2 α j = 0 s j = 0 s j = 0 x j (17.2) (17.14c) p j c j p j = 0 j c j p T k S j = k 0 (17.13)(e) (f) β n 1 = 0, p n = 0 (17.13) α k β k (17.10) (17.8) p T k s k = (s k + β k 1 p k 1 ) T s k = s T k s k (17.15) (17.6) (17.13)(b) α k (17.4) Sp k = 1 S(x k+1 x k ) = 1 (s k s k+1 ) α k α k s T k+1sp k = 1 s T α k+1(s k s k+1 ) k = 1 s k+1 2 = p T s k+1 2 k Sp α k k s k 2 (17.12) (17.13)(e) β k s k+1 β k+1 0 α k S p k = 0 0 p k = 0 (17.15) s k = 0 s k+1 (CGNE) A 0 = j c j p T k Sp j = c k p T k Sp k Ax = b (17.16) c k 0 p k n x n = x 0 + α 0 p 0 + α 1 p 1 + + α n 1 p n 1 n n α k (17.2) s k n n s n = 0 x n (17.2) CG A A T Ax = A T b (17.17) (17.16) (17.17) (17.2) S = A T A c = A T b (17.18) (17.3) F (x) = 1 2 Ax b 2 1 2 b 2 (17.19) A Ax b 2 = min (17.20)

17. 4 (17.13) (17.18) (17.13) A T A A p T k Sp k = p T k A T Ap k = Ap k 2 A T A CGNE (a) r 0 = b Ax 0 p 0 = s 0 = A T r 0 (b) q k = Ap k (c) α k = s k 2 q k 2 (d) x k+1 = x k + α k p k (17.21) (e) r k+1 = r k α k q k (f) s k+1 = A T r k+1 (g) β k = s k+1 2 s k 2 (h) p k+1 = s k+1 + β k p k (b) r k (17.16) r k = b Ax k (17.22) α k r k+1 (17.21)(e) r k+1 2 = r k 2 2α k q T k r k + α 2 k q k 2 r k+1 2 α k α k 0 α k = qt k r k q k 2 q T k r k = (Ap k ) T r k = p T k A T r k = p T k s k α k (17.6) α k CGNE r k+1 r k+1 2 = r k 2 α 2 k q k 2 (17.23) 1 p k s k (17.14a) (17.14b) (17.14a) p T j s k = q T j r k = 0 j < k (17.14c) p T j A T Ap k = q T j q k = 0 j k s k q k CGNE q T j r k = 0 j < k (17.24a) s T j s k = 0 j k (17.24b) q T j q k = 0 j k (17.24c) CGNE s k+1 α k q k 2 0 = q T k r k = p T k s k = (s k + β k 1 p k 1 ) T s k (17.14a) 0 q k 0 s k 0 s k+1 = 0 0 CGNE A A n m (17.16) (17.20) x n > m n m x k, p k, s k : m r k, q k : n

17. 5 s k q k n > m s m 0 s k s m = A T r m = A T (b Ax m ) = 0 x m A x 0 ( ) A (17.21) normal equation CGNE A A CGNE x 0 p k A T p k R(A T ) p 0 = A T r 0 p 0 p 1 (17.21) q 0 = Ap 0 r 1 = r 0 α 0 Ap 0 s 1 = p 0 α 0 A T Ap 0 p 1 = (1 + β 0 )p 0 α 0 A T Ap 0 p 1 R(A T ) p k (17.21) p k+1 = (1 + β k )p k α k A T Ap k β k 1 p k 1 p 0, p 1 p k p k = {p 0, Bp 0, B 2 p 0,, B k p 0 } R(A T ) B = A T A x k k 1 x k = x 0 + α j p j j=0 CGNE x 0 R(A T ) x k R(A T ) 3 CGNE x 0 R(A T ) x = A b (17.25) A A A T (b Ax) = 0 x 0 R(A T ) x k n y 0 x 0 = A T y 0 x 0 A x 0 0 x k R(A T ) (17.16) b x 0 b = b + b b N(A T ) b R(A) x 0 = x 0 + x 0 (17.26) x 0 N(A) x 0 R(A T ) A T b = 0 Ax 0 = 0 x k+1 [ k ] x k+1 = x 0 + x 0 + α j p j j=0 N(A) R(A T ) N(A) x k+1 [ r k+1 = b + (b Ax 0) k ] α j q j j=0 N(A T ) q j R(A) R(A) b

17. 6 (CGMN) (17.16) AA T z = b x = A T z (17.27) x 2 ( 10, (10.52) ) CGS (17.13) B = AA T z k x k CGMN (a) p 0 = r 0 = b Ax 0 (b) q k = A T p k (c) α k = r k 2 q k 2 (d) x k+1 = x k + α k q k (17.28) (e) r k+1 = r k α k Aq k (f) β k = r k+1 2 r k 2 (g) p k+1 = r k+1 + β k p k (b) r k (17.22) p k q k (17.16) b b R(A) ( 11) h Ah = b (17.29) x k e k = h x k (17.30) (17.28)(d) e k+1 = e k α k q k (17.31) e k+1 e k+1 2 = e k 2 2α k q T k e k + α 2 k q k 2 q k q T k e k = p T k Ae k = p T k (Ah Ax k ) = p T k r k e k+1 2 = e k 2 2α k p T k r k + α 2 k q k 2 α k = pt k r k q k 2 (17.32) α k e k+1 2 = e k 2 α 2 k q k 2 (17.33) 1 b A b β k q T k q k+1 = 0 (17.34) x k (17.28) q k+1 = A T r k+1 + β k q k 0 = q T k+1q k = r T k+1aq k + β k q k 2 β k = rt k+1 Aq k q k 2 (17.35) p j, q j, r j CGNE p T j r k = 0 j < k (17.36a) r T j r k = 0 j k (17.36b) q T j q k = 0 j k (17.36c) (17.32) α k (17.35) β k (17.28) (17.28)(g) (17.28)(a) p T k r k = r k 2 + β k 1 p T k 1r k = r k 2 (17.32) (13.28)(c) (17.28) (e) (17.28)(b) r k+1 2 = r T k+1(r k α k Aq k ) = α k r T k+1aq k

17. 7 α k (17.35) β k (17.28)(f) CGMN r k = 0 q k = 0 r k = 0 (17.16) r j n n 0 CGNE r j r k 0 r j 0 q k = 0 CGNE r k = 0 b R(A) q k = 0 r k = 0 (17.31) q T k e k+1 = q T k e k α k q k 2 (17.28)(b) q T k e k = p T k Ae k = p T k r k (17.30) Ae k = r k q T k e k+1 = p T k r k α k q k 2 = 0 (17.32) e T k q k = e T k (A T r k + β k 1 q k 1 ) = r k 2 + β k 1 e T k q k 1 = r k 2 q k = 0 r k = 0 (17.28)(b) q k R(A T ) x 0 R(A T ) x k R(A T ) b R(A) x 0 R(A T ) 4 CGMN x = A b Ax = b ( 10 ) minimum norm CGMN z b R(A) x b R(A) (17.32) h = A b + x 0 Ax 0 = 0 e k x k R(A T ) e k 2 = A b + x 0 x k 2 = A b x k 2 + x 0 2 (17.33) 1 A b x k 2 1 (17.2) (17.16) S (17.2) C SC 1 x = b x = Cx C 1 Sx = C 1 b C κ(s) κ(sc 1 ) κ(c 1 S) C = S (17.2) C S C C S (17.17) Sz = c z = Cx (17.37) S = C T A T AC 1 c = C T A T b C T (C 1 ) T C A T A A T AC 1 I

17. 8 C Ct = p (17.38) A n m C m m (17.37) S CGS(17.13) z x GCPCNE (a) r 0 = b Ax 0 p 0 = s 0 = C T A T r 0 (b) q k = AC 1 p k (c) α k = s k 2 q k 2 (d) x k+1 = x k + α k C 1 p k (17.39) (e) r k+1 = r k α k q k (f) s k+1 = C T A T r k+1 (g) β k = s k+1 2 s k 2 (h) p k+1 = s k + β k p k PC (preconditioned) C 1 C 1 p (17.38) CGMN Björck, Å. and T. Elfving (1979) : Accelerated projection methods for computing pseudoinverse solutions of systems of linear equations, BIT, 19, 145-163. R(A) x y = Ax y R(A) A (range) N(A) Az = 0 z N(A) A (null space) R(A T ) N(A) A n m m R(A T ) N(A) A A = UΛV T V m m m x V v i 0 v i Av i = 0 N(A) 0 U u i λ i A T u i = λ iv i λ i > 0 R(A T ) 1 x = A b b Ax x x 11 2 (i) x R(A T ) (ii) (b Ax) R(A) (i) x N(A) (ii) (b Ax) N(A T ) x = A b (ii) x A T (b Ax) = 0 x x Ax = 0 x (i) x N(A) 3 A T Ax = A T b x R(A T ) 4 x = A b Ax = b b R(A) x R(A T ) x = A b