Bài tập lớn môn Phương pháp tính PDF

Preview

Summary

ĐẠI HỌC QUỐC GIA THÀNH PHỐ HỒ CHÍ MINHTRƯỜNG ĐẠI HỌC BÁCH KHOABÁO CÁO BÀI TẬP LỚN MÔN HỌC####### PHƯƠNG PHÁP TÍNHCHỦ ĐỀ 3LỚP: L07 – NHÓM: 17 – HKGVHD: TS. NGUYỄN ĐÌNH DƯƠNGTP. HỒ CHÍ MINH, NĂM 2021SINH VIÊN THỰC HIỆNSTT MSSV Họ Tên Điểm BTL1 2012813 Đinh Thái Nhật Duy2 2013184 Đoàn Minh Hiể...

Description

ĐẠI HỌC QUỐC GIA THÀNH PHỐ HỒ CHÍ MINH TRƯỜNG ĐẠI HỌC BÁCH KHOA



BÁO CÁO BÀI TẬP LỚN MÔN HỌC

PHƯƠNG PHÁP TÍNH

CHỦ ĐỀ 3 LỚP: L07 – NHÓM: 17 – HK211 GVHD: TS. NGUYỄN ĐÌNH DƯƠNG

SINH VIÊN THỰC HIỆN STT

MSSV

1

2012813

Đinh Thái Nhật

Duy

2

2013184

Đoàn Minh

Hiển

3

2013322

Nguyễn Thái

Huy

4

2013868

Lý Gia

Nghĩa

5

1814485

Nguyễn Đức

Trọng

Họ

Tên

TP. HỒ CHÍ MINH, NĂM 2021

Điểm BTL

BẢNG PHÂN CÔNG NHIỆM VỤ VÀ THEO DÕI TIẾN ĐỘ LÀM VIỆC NHÓM

Nhiệm vụ được phân công

STT

MSSV

1

2012813

Đinh Thái Nhật Duy

Giải bài 3 + thuyết trình bài 3

2

2013184

Đoàn Minh Hiển

Bài 2 + code bài 2 + slides bài 2 + thuyết trình bài 2

3

2013322

Nguyễn Thái Huy

Bài 1 + code bài 1 + code bài 3 + slides bài 3 + tổng hợp + thuyết trình bài 1

4

2013868

Lý Gia Nghĩa

5

1814485

Nguyễn Đức Trọng

Họ và tên

– Lý thuyết bài 3 + slides bài 3

Điểm BTL

Ký tên

MỤC LỤC PROBLEM 1 ............................................................................................................................. 1 1.1

1.2

1.3

1.4

Lý thuyết ...................................................................................................................... 1 1.1.1

Phương trình phi tuyến ..................................................................................... 1

1.1.2

Phương trình vi phân ........................................................................................ 3

Nội dung và kết quả có được ....................................................................................... 5 1.2.1

Nội dung ............................................................................................................ 5

1.2.2

Kết quả .............................................................................................................. 6

Bài toán mở rộng ....................................................................................................... 12 1.3.1

Nội dung .......................................................................................................... 12

1.3.2

Kết quả ............................................................................................................ 13

Kết luận ...................................................................................................................... 15

PROBLEM 2 ........................................................................................................................... 16 2.1

Lý thuyết .................................................................................................................... 16

2.2

Nội dung và kết quả có được ..................................................................................... 19 2.2.1

Nội dung .......................................................................................................... 19

2.2.2

Kết quả ............................................................................................................ 20

2.3. Bài toán mở rộng ....................................................................................................... 23 2.3.1.

Nội dung .......................................................................................................... 23

2.3.2.

Kết quả ............................................................................................................ 23

2.4. Kết luận ...................................................................................................................... 24 PROBLEM 3 ........................................................................................................................... 25 3.1

3.2

3.3

3.4

Lý thuyết .................................................................................................................... 25 3.1.1

Phương pháp Euler cải tiến cho hệ phương trình vi phân ............................. 25

3.1.2

Spline bậc ba tự nhiên..................................................................................... 25

Nội dung và kết quả có được ..................................................................................... 27 3.2.1

Nội dung .......................................................................................................... 27

3.2.2

Kết quả ............................................................................................................ 28

Bài toán mở rộng ....................................................................................................... 37 3.3.1

Nội dung .......................................................................................................... 37

3.3.2

Kết quả ............................................................................................................ 37

Kết luận ...................................................................................................................... 38

TÀI LIỆU THAM KHẢO...................................................................................................... 40

PROBLEM 1 1.1 Lý thuyết Cho f ( x) là hàm số liên tục trên đoạn  a, b thỏa mãn f ( a) f ( b )  0 . Khi đó, tồn tại một giá trị p nằm giữa a và b sao cho f ( p ) = 0 . Giả sử giá trị p là duy nhất, ta gọi  a, b là khoảng cách ly nghiệm của phương trình f ( x ) = 0 . Gọi x * là nghiệm gần đúng của nghiệm đúng p trên  a, b thì ta có công thức đánh giá sai số tổng quátEquation Chapter 1 Section 1 x* − p 

f ´( x *) , m

m = min f ´( x)  0 . x  a , b 

(1.1)

1.1.1 Phương trình phi tuyến 1.1.1.1 Phương pháp chia đôi Nội dung: Kiểm tra giá trị của hàm tại điểm giữa khoảng c =

a +b . Nếu f ( c ) f ( a )  0 thì 2

chắc chắn rằng nghiệm của phương trình sẽ nằm trong  a, c có độ dài bằng một nửa đoạn  a, b ban đầu. Ngược lại, nếu f ( c ) f ( b )  0 thì ta có điều tương tự đối với đoạn c, b  . Lặp lại các bước này sẽ giúp định vị ngày càng chính xác hơn vị trí nghiệm của phương trình. Đây gọi là phương pháp chia đôi (the bisection method). Công thức đánh giá sai số: Giả sử sau n lần chia đôi, ta tìm được đoạn  an , bn  có độ dài nghiệm gần đúng là điểm giữa khoảng xn =

b−a . Chọn 2n

an + bn vì nó gần với giá trị nghiệm chính 2

xác nhất, khi đó ta có công thức đánh giá sai số:

xn − p 

1

b−a . 2n +1

(1.2)

1.1.1.2 Phương pháp cát tuyến Ở phương pháp Newton (Newton’s method hay Newton-Raphson method), từ điểm có hoành độ x n− 1 trên đồ thị của đường cong y = f ( x ) , ta kẻ tiếp tuyến với đường cong. Hoành độ giao điểm của tiếp tuyến với trục hoành sẽ là xn . Ta dễ dàng viết phương trình tiếp tuyến:

y = f ´( x n−1 )( x − x n−1 ) + f ( x n−1 ) . Cho y = 0, x = xn ta thu được công thức xác định xn , cũng chính là công thức lặp của phương pháp Newton:

xn = xn − 1 −

f ( x n −1 )

f ´ ( x n− 1 )

, n = 1, 2, 3,...

(1.3)

Đây là một trong những công thức nổi tiếng và được sử dụng rộng rãi nhất trong việc giải quyết vấn đề tìm nghiệm. Tuy nhiên, nó có một điểm yếu, đó là phải biết được giá trị đạo hàm của f ở mỗi lần xấp xỉ. Đôi khi, việc tìm f ´( x ) sẽ rất khó khăn và không thuận tiện do cần nhiều kỹ thuật tính toán. Với những trường hợp này, ta có thể tính gần đúng đạo hàm bằng công thức

f ´( x n ) 

f ( xn ) − f ( xn − 1 ) x n − x n− 1

,

khi đó, phương trình (1.3) có thể được viết lại thành công thức lặp như sau:

xn = xn − 1 −

f ( xn −1 )( xn −1 − xn − 2 ) f ( x n −1 ) − f ( x n − 2 )

, n = 2, 3, 4,...

(1.4)

Kỹ thuật này được gọi là phương pháp cát tuyến (the secant method) (xem Hình 1.1). Xuất phát từ hai giá trị gần đúng ban đầu p0 , p1 , giá trị p2 là giao điểm của trục x với đường thẳng đi qua hai điểm ( p 0, f ( p 0 ) ) và ( p 1, f ( p 1) ) . Giá trị p3 là giao điểm của trục x với đường thẳng đi qua hai điểm ( p 1, f ( p 1) ) và ( p 2, f ( p 2 ) ) , tương tự với các giá trị xấp xỉ tiếp theo.

2

Hình 1.1 Minh họa phương pháp cát tuyến Khi muốn đánh giá sai số nghiệm gần đúng đã tìm được, ta dùng công thức ước lượng sai số tổng quát (1.1). 1.1.2 Phương trình vi phân Cho bài toán Cauchy:

 y ´( x ) = f ( x, y ( x )) ,    y( x ) = y 0 0  

x  x0 ,

(1.5)

với y = y ( x ) là hàm cần tìm, khả vi trên đoạn  a, b , y0 là giá trị ban đầu cho trước của y ( x ) tại x = x 0 , f ( x, y) là hàm hai biến liên tục với các đạo hàm riêng của nó. Để tìm nghiệm gần đúng của bài toán (1.5), ta chia đoạn  a, b thành n đoạn nhỏ bằng nhau với bước chia h =

b−a . Khi đó, các điểm chia là x0 = a, xk = x0 + kh, n

k = 0,1,2,..., n, xn = b . 1.1.2.1 Phương pháp Euler cải tiến Giả sử y ( x ) là nghiệm duy nhất của bài toán (1.5), có đạo hàm cấp hai liên tục trên  a, b . Ta xây dựng đượ c công thức Euler (Euler’s method) như sau:

y ( xk +1 )  yk +1 = yk + hf ( xk , yk ), k = 0,1, 2,..., n −1 3

(1.6)

f ( xk , yk ) + f ( xk +1 , yk +1 ) ta được 2

Nếu thay f ( x k, y k ) trong công thức (1.6) bởi công thức Euler cải tiến (modified Euler’s method):

y ( xk +1 )  yk +1 = yk + h

f ( x k , y k ) + f ( x k+1 , y k+1 ) , k = 0,1,2,..., n −1 . 2

Để đơn giản, ta thay yk +1 ở vế phải bởi yk −1 + hf ( xk −1 , yk−1 ) , khi đó:

y ( xk +1 )  yk +1 = yk + h

f ( xk , y k ) + f ( xk+ 1, yk + hf ( xk , yk ) ) 2

, k = 0,1,2,..., n − 1 .

Viết lại công thức Euler cải tiến:

 k  K = hf ( x , y )  1 k k   k  K = hf ( x + h , y + K k ) 1 k k  2   k k   y ( x )  y = y + K1 + K2  k +1 k +1 k 2 

(1.7)

1.1.2.2 Phương pháp Runge-Kutta bậc bốn

  y = y (x + h )  y  k +1 k k    K k = hf ( x , y )  1 k k  k  K 2 = hf ( xk +  2h , y k   k  K = hf x +  h , y (k 3 k  3   M   k  K n = hf ( x k + nh , y k 

n

+ Aj K jk j =1

+ 21K 1k ) + 31K + 32K k 1

(1.8) k 2

)

+ n1K1k + n2 K 2k + + n, n −1K nk−1 )

Ta cần xác định các hệ số A1, A2, K, An;  2,  3, K,  n;  21,  31, K, n,n −1 . Khai triển Taylor nghiệm y ( x ) tại xk đến bậc m rồi thay x = xk +1 , ta được:

y ( xk +1 ) = y ( xk + h )  y ( xk ) + y ´( xk ) h + y ˝ ( xk 3

4

)

h2 2!

h h hm ( 4) (m) + y ( x k ) + y ( x k ) + K + y ( xk ) 3! 4! m! ( 3)

4

Trong trường hợp m = n = 4 , ta có công thức Runge-Kutta bậc bốn (FourthOrder Runge-Kutta method – RK4):

  y k + = y ( xk + h )  y k + 1 ( K k + 2 K k + 2 K k + K k ) 2 3 4  1 6 1   k  K1 = hf ( xk , yk )    h Kk  k  K 2 = hf  xk + , y k + 1    2 2     k   k  K = hf  x + h , y + K 2   3  k 2 k 2      k k  K 4 = hf ( xk + h , yk + K 3 ) 

(1.9)

1.2 Nội dung và kết quả có được 1.2.1 Nội dung A bungee jumper jumps from a mountain with the downward vertical velocity v described by the mathematical model:

c dv = g − d v2 (see the picture), where m is the dt m

mass of jumper and cd is called drag coefficient. a) Suppose that the jumper is initial at rest, find analytically the expression of v. b) Let g = 9.8 (m/s 2), m = 68.1 (kg), c d = 0.25 (kg/m) and the jumper is initial at rest, establish the table to compute the velocity of the jumper for the first 10 seconds with step size h = 1 (s) by using modified Euler’s and Runge-Kutta method. Compare the results to the exact values found in a). c) Using the result of a) and the bisection method, the secant method to determine the drag coefficient for a jumper with the weight of 95 (kg) and the velocity

v = 46 (m/s) after 10 seconds of fall until the relative error is less than 5% (Guess the isolated interval containing root). Tạm dịch: Một vận động viên bungee nhảy từ một ngọn núi với vận tốc v thẳng đứng hướng xuống (xem hình), được mô tả bởi mô hình toán học 5

dv c = g − d v2 dt m

(1.10), trong đó,

m là khối lượng của vận động viên, cd là hệ số cản.

a) Giả sử ban đầu vận động viên đang ở trạng thái nghỉ, tìm biểu thức mô tả v . b) Biết g = 9.8 m/s2, m = 68.1 (kg), cd = 0.25 (kg/m) và vận động viên bắt đầu nhảy từ trạng thái nghỉ. Lập bảng để tính vận tốc của vận động viên trong 10 giây đầu tiên kể từ lúc rời ngọn núi với bước chia h = 1(s) bằng phương pháp Euler cải tiến và phương pháp Runge-Kutta. So sánh kết quả với giá trị chính xác tìm được từ câu a). c) Sử dụng kết quả câu a) cùng các phương pháp chia đôi, phương pháp cát tuyến để tìm hệ số cản đối với vận động viên có khối lượng 95 (kg), vận tốc v = 46 (m/s) sau 10 giây kể từ lúc bắt đầu rơi đến khi sai số tương đối nhỏ hơn 5% (dự đoán khoảng cách ly nghiệm). 1.2.2 Kết quả a) Nhân hai vế của phương trình (1.10) với

m , ta được cd

m dv gm 2 = −v . cd d t cd Đặt a =

gm , phương trình trở thành cd m dv = a2 − v2 cd d t



6

dv c = d dt 2 a −v m 2

 

dv c =  d dt 2 a −v m  v c 1 artanh   = d t + C a  a m



2

Do lúc đầu, vận động viên ở trạng thái nghỉ nên v ( 0 ) = 0 , suy ra C = 0 . 

 c  v = a tanh  a d t  m 



v=

 gc  gm  d  t tanh    cd m  

(1.11)

b) Giá trị vận tốc chính xác vexact được tính bằng biểu thức (1.11), với g = 9.8,

cd = 0.25, m = 68.1, kết quả được thể hiện ở cột 2 của Bảng 1.1. Với phương pháp Euler cải tiến, từ công thức (1.7): Khi t = 0:

  0.25 K10 = 1   9.8 −  02  = 9.8 68.1    2 0.25 K 20 = 1   9.8 −  ( 0 + 9.8 )   9.4474303 68.1   + 9.8 9.4474303 vm.Euler 1  0 +  9.6237151 2 vm.Euler 1 − vexact 1  0.0604286 Khi t = 1 :

  0.25 K11 = 1   9.8 −  9.62371512   9.4600004 68.1     0.25 2 1  (9.6237151 + 9.4600004 )   8.4630389 K 2 = 1   9.8 − 68.1   9.4600004 + 8.4630389 vm.Euler 2  9.6237151 +  18.5852348, 2 vm.Euler 2 − vexact 2  0.1257200

7

Thao tác tương tự đối với các giá trị t tiếp theo, ta thu được kết quả ở cột 3 và cột 5 của Bảng 1.1. Trong đó, v m.Euler 10  49.2271057 , vm.Euler 10 − v  0.1647634 . Với phương pháp Runge-Kutta bậc bốn (RK4), áp dụng công thức (1.9): Với t = 0:

  0.25 K10 = 1   9.8 −  02  = 9.8 68.1   2  0.25  9.8    0 K 2 = 1   9.8 −   9.7118576 0+  68.1  2    2  0.25  9.7118576    0 K 30 = 1   9.8 − +     9.7134360  68.1 2      0.25 2 K 40 = 1   9.8 − ( 0 + 9.7134360)   9.45 36313 68.1   1 vRK41  0 + (9.8 + 2  9.7118576 + 2  9.7134360 +9.4536313)  9.6840364, 6 vRK4 1 − vexact 1  0.0001073 Với t = 1 :

  0.25 K11 = 1   9.8 −  9.68403642  = 9.4557248 68.1   2  0.25  9.4557248    1 K 2 = 1   9.8 −  9.6840364 +    9.0375080  68.1 2     2  0.25  9.0375080    1 K 3 = 1   9.8 −  9.6840364 +    9.0594741  68.1 2      0.25 2 K 14 = 1   9.8 − (0 + 9.0594741)   8.5102820   68.1 1 vRK4 2  0 + (9.4557248 + 2  9.0375080 + 2  9.0594741 + 8.5102820 ) = 18.7106982 6 vRK4 2 − vexact 2  0.0002565

8

Thao tác tương tự với các giá trị t còn lại, ta thu được kết quả ở cột 4 và cột 6 của Bảng 1.1. Trong đó, vRK410  49.3909717 và vm.Euler 10 − vexact 10  0.0008974 . Điều đó cho thấy rằng việc sử dụng phương pháp Runge-Kutta bậc bốn mặc dù phức tạp hơn nhưng cho giá trị xấp xỉ tốt hơn nhiều so với phương pháp Euler cải tiến. Bảng 1.1 So sánh giá trị vận tốc được tính theo phương pháp Euler cải tiến và phương pháp Runge-Kutta bậc bốn với giá trị chính xác t

vexact

vm.Euler

vRK4

vm.Euler − v

vRK4 − v

0

0

0

0

0

0

1

9.6841437

9.6237151

9.6840364

0.0604286

0.0001073

2

18.7109547

18.5852348

18.7106982

0.1257200

0.0002565

3

26.5902279

26.4014722

26.5897509

0.1887557

0.0004770

4

33.0831500

32.8445302

33.0824094

0.2386199

0.0007406

5

38.1845782

37.9175859

38.1835979

0.2669923

0.0009803

6

42.0446493

41.7725478

42.0435124

0.2721015

0.0011369

7

44.8830425

44.6250828

44.8818546

0.2579597

0.0011879

8

46.9265585

46.6953182

46.9254134

0.2312403

0.0011451

9

48.3755332

48.1770315

48.3744953

0.1985017

0.0010379

10

49.3918691

49.2271057

49.3909717

0.1647634

0.0008974

c) Thay các giá trị: g = 9.8, m = 95, v = 46, t = 10 vào phương trình (1.11) ta được:

 9.8 c  9.8  95   d tanh  f ( cd ) =  10 − 46 = 0   95 cd   Kiểm tra điều kiện phương trình (1.12) trong đoạn 0.3, 0.5  :

 9.8  0.3  9.8  95   tanh  f (0.3 ) =  10  − 46  6.4995   0.3 95    9.8  0.5  9.8  95   tanh   10  − 46  −3.7585   0.5 95   f ( 0.3) f ( 0.5)  0.

f (0.5 ) = 

9

(1.12)

Ngoài ra, dùng Maple để kiểm tra, ta thấy: f ´ (cd )  0 cd  0.3, 0.5  , do đó,

0.3, 0.5  là khoảng cách ly nghiệm của phương trình (1.12). Phương pháp chia đôi:

0.3 + 0.5 = 0.4 2 f ( 0.4 )  0.6123466  0  a 1 = c d 0 = 0.4;

cd0 =

0.4 + 0.5 = 0.45 2 f ( 0.45 )  −1.7218725  0  b2 = c d 1 = 0.45.

cd1 =

Sai số tương đối của cd1 là 1 =

cd1 − cd 0 cd 1

 100% =

0.45 − 0.4  11.1111111%. 0.45

Quá trình tính toán tiếp diễn với các bước tương tự được thể hiện trong Bảng 1.2. Bảng 1.2 Hệ số cản tính theo phương pháp chia đôi

n (%)

( )

n

an

bn

cdn

0

0.3

0.5

0.4

1

0.4

0.5

0.45

11.1111111

–1.7218725

2

0.4

0.45

0.425

5.8823529

–0.5960967