Skip to content Skip to navigation Skip to collection information

OpenStax-CNX

You are here: Home » Content » Giáo trình giải tích mạng điện » Giải phương trình vi phân bằng phương pháp số

Navigation

Lenses

What is a lens?

Definition of a lens

Lenses

A lens is a custom view of the content in the repository. You can think of it as a fancy kind of list that will let you see content through the eyes of organizations and people you trust.

What is in a lens?

Lens makers point to materials (modules and collections), creating a guide that includes their own comments and descriptive tags about the content.

Who can create a lens?

Any individual member, a community, or a respected organization.

What are tags? tag icon

Tags are descriptors added by lens makers to help label content, attaching a vocabulary that is meaningful in the context of the lens.

This content is ...

Affiliated with (What does "Affiliated with" mean?)

This content is either by members of the organizations listed or about topics related to the organizations listed. Click each link to see a list of all content affiliated with the organization.
  • VOCW

    This collection is included inLens: Vietnam OpenCourseWare's Lens
    By: Vietnam OpenCourseWare

    Click the "VOCW" link to see all content affiliated with them.

Recently Viewed

This feature requires Javascript to be enabled.
 

Giải phương trình vi phân bằng phương pháp số

Module by: PGS. TS. Lê Kim Hùng. E-mail the author

2.1. GIỚI THIỆU.

Nhiều hệ thống vật lý phức tạp được biểu diễn bởi phương trình vi phân nó không có thể giải chính xác bằng giải tích. Trong kỹ thuật, người ta thường sử dụng các giá trị thu được bằng việc giải gần đúng của các hệ phương trình vi phân bởi phương pháp số hóa. Theo cách đó, lời giải của phương trình vi phân đúng là một giai đoạn quan trọng trong giải tích số.

Trong trường hợp tổng quát, thứ tự của việc làm tích phân số là quá trình từng bước chính xác chuổi giá trị cho mỗi biến phụ thuộc tương ứng với một giá trị của biến độc lập. Thường thủ tục là chọn giá trị của biến độc lập trong một khoảng cố định. Độ chính xác cho lời giải bởi tích phân số phụ thuộc cả hai phương pháp chọn và kích thước của khoảng giá trị. Một số phương pháp thường xuyên dùng được trình bày trong các mục sau đây.

2.2. GIẢI PHƯƠNG TRÌNH VI PHÂN BẰNG PHƯƠNG PHÁP SỐ.

2.2.1 Phương pháp Euler:

Cho phương trình vi phân bậc nhất.

dydx=f(x,y)dydx=f(x,y) size 12{ { { ital "dy"} over { ital "dx"} } =f \( x,y \) } {} (2.1)

yxyxy = g(x,c)y0x0Hình 2.1: Đồ thị của hàm số từ bài giải phương trình vi phân0

Khi x là biến độc lập và y là biến phụ thuộc, nghiệm phương trình (2.1) sẽ có dạng:

y = g(x,c) (2.2)

Với c là hằng số đã được xác định từ lý thuyết trong điều kiện ban đầu. Đường cong miêu tả phương trình (2.2) được trình bày trong hình (2.1). Từ chỗ tiếp xúc với đường cong, đoạn ngắn có thể giả sử là một đoạn thẳng. Theo cách đó, tại mỗi điểm riêng biệt (x0,y0) trên đường cong, ta có:

Δy dy dx 0 Δx Δy dy dx 0 Δx size 12{Δy approx { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } Δx} {}

Với dydx0dydx0 size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } } {}là độ dốc của đường cong tại điểm (x0,y0). Vì thế, ứng với giá trị ban đầu x0 và y0, giá trị mới của y có thể thu được từ lý thuyết là x:

y1=y0+Δyy1=y0+Δy size 12{y rSub { size 8{1} } =y rSub { size 8{0} } +Δy} {} hay y1=y0+dydx0hy1=y0+dydx0h size 12{y rSub { size 8{1} } =y rSub { size 8{0} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } h} {} (đặt h = x)

Khi y là số gia của y tương ứng với một số gia của x. Tương tự, giá trị thứ hai của y có thể xác định như sau.

y 2 = y 1 + dy dx 1 h y 2 = y 1 + dy dx 1 h size 12{y rSub { size 8{2} } =y rSub { size 8{1} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{1} } h} {}

xx0x1x2x3y0y1y2y3hhhy= g(x,c)Hình 2.2 : Đồ thị của lời giải xấp xỉcho phương trình vi phân bằng phương pháp Euler0y

Khi dydx1=f(x1,y1)dydx1=f(x1,y1) size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{1} } =f \( x rSub { size 8{1} } ,y rSub { size 8{1} } \) } {}

Quá trình có thể tính tiếp tục, ta được:

y 3 = y 2 + dy dx 2 h y 3 = y 2 + dy dx 2 h size 12{y rSub { size 8{3} } =y rSub { size 8{2} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{2} } h} {}

y 4 = y 3 + dy dx 3 h y 4 = y 3 + dy dx 3 h size 12{y rSub { size 8{4} } =y rSub { size 8{3} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{3} } h} {}

...........................

Bảng giá trị x và y cung cấp cho toàn bộ bài giải phương trình (2.1). Minh họa phương pháp như hình 2.2.

2.2.2. Phương pháp biến đổi Euler.

Trong khi ứng dụng phương pháp Euler, giá trị dy/dx của khoảng giả thiết tính toán bắt đầu vượt ra ngoài khoảng cho phép. Sự thay thế đó có thể thu được bằng cách tính toán giá trị mới của y cho x1 như trước.

x1 = x0 + h

y 1 ( 0 ) = y 0 + dy dx 0 h y 1 ( 0 ) = y 0 + dy dx 0 h size 12{y rSub { size 8{1} } rSup { size 8{ \( 0 \) } } =y rSub { size 8{0} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } h} {}

Dùng giá trị mới x1 và y1(0) thay vào phương trình (2.1) để tính toán gần đúng giá trị của dydx1dydx1 size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{1} } } {}tại cuối khoảng.

dy dx 1 ( 0 ) = f ( x 1 , y 1 ( 0 ) ) dy dx 1 ( 0 ) = f ( x 1 , y 1 ( 0 ) ) size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{1} } rSup { size 8{ \( 0 \) } } =f \( x rSub { size 8{1} } ,y rSub { size 8{1} } rSup { size 8{ \( 0 \) } } \) } {}

Sau đó tận dụng giá trị y1(1) có thể tìm thấy bởi dùng trung bình của dydx0dydx0 size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } } {}dydx1(0)dydx1(0) size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{1} } rSup { size 8{ \( 0 \) } } } {}như sau:

y 1 ( 1 ) = y 0 + dy dx 0 + dy dx 1 ( 0 ) 2 h y 1 ( 1 ) = y 0 + dy dx 0 + dy dx 1 ( 0 ) 2 h size 12{y rSub { size 8{1} } rSup { size 8{ \( 1 \) } } =y rSub { size 8{0} } + left ( { { { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{1} } rSup { size 8{ \( 0 \) } } } over {2} } right )h} {}

Dùng x1 và y1(1), giá trị xấp xỉ thứ ba y1(2) có thể thu được bởi quá trình tương tự như sau: y1(2)=y0+dydx0+dydx1(1)2hy1(2)=y0+dydx0+dydx1(1)2h size 12{y rSub { size 8{1} } rSup { size 8{ \( 2 \) } } =y rSub { size 8{0} } + left ( { { { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{1} } rSup { size 8{ \( 1 \) } } } over {2} } right )h} {}

Ta được:

y 1 ( 3 ) = y 0 + dy dx 0 + dy dx 1 ( 2 ) 2 h y 1 ( 3 ) = y 0 + dy dx 0 + dy dx 1 ( 2 ) 2 h size 12{y rSub { size 8{1} } rSup { size 8{ \( 3 \) } } =y rSub { size 8{0} } + left ( { { { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{1} } rSup { size 8{ \( 2 \) } } } over {2} } right )h} {}

Quá trình có thể tính tiếp tục cho đến khi hai số liền nhau ước lượng cho y là ngang bằng nằm trong phạm vi mong muốn. Quá trình hoàn toàn lặp lại thu được giá trị y2. Kết quả thu được có sự chính xác cao hơn từ sự biến đổi của phương pháp Euler được minh họa trong hình 2.3.

Hình 1
Hình 1 (.wmf)
y = g(x,c)yxx0x1hy0
Hình 2
Hình 2 (.wmf)
Hình 2.3 : Đồ thị của lời giải xấp xỉ cho phương trình vi phân bằng phương pháp biến đổi Euler.0y1y2dy (0)dx 1

Phương pháp Euler có thể ứng dụng để giải hệ phương trình vi phân cùng lúc. Cho hai phương trình:

dy dx = f 1 ( x ,y,z ) dz dx = f 2 ( x ,y,z ) dy dx = f 1 ( x ,y,z ) dz dx = f 2 ( x ,y,z ) alignl { stack { size 12{ { { ital "dy"} over { ital "dx"} } =f rSub { size 8{1} } \( x",y,z" \) } {} # { { ital "dz"} over { ital "dx"} } =f rSub { size 8{2} } \( x",y,z" \) {} } } {}

Với giá trị ban đầu x0, y0 và z0 giá trị mới y1 sẽ là:

y 1 = y 0 + dz dx 0 h y 1 = y 0 + dz dx 0 h size 12{y rSub { size 8{1} } =y rSub { size 8{0} } + { { ital "dz"} over { ital "dx"} } \rline rSub { size 8{0} } h} {}

Với: dydx0=f1(x0,y0,z0)dydx0=f1(x0,y0,z0) size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } =f rSub { size 8{1} } \( x rSub { size 8{0} } ",y" rSub { size 8{0} } ",z" rSub { size 8{0} } \) } {}

Tương tự.

z 1 = z 0 + dz dx 0 h z 1 = z 0 + dz dx 0 h size 12{z rSub { size 8{1} } =z rSub { size 8{0} } + { { ital "dz"} over { ital "dx"} } \rline rSub { size 8{0} } h} {}

Với: dzdx0=f2(x0,y0,z0)dzdx0=f2(x0,y0,z0) size 12{ { { ital "dz"} over { ital "dx"} } \rline rSub { size 8{0} } =f rSub { size 8{2} } \( x rSub { size 8{0} } ,y rSub { size 8{0} } ,z rSub { size 8{0} } \) } {}

Cho số gia tiếp theo, giá trị x1 = x0 + h, y1 và z1 dùng để xác định y2 và z2. Trong phương pháp biến đổi Euler y1 và z1 dùng để xác định giá trị đạo hàm tại x1 cho đánh giá gần đúng cấp hai y1(1) và z1(1).

2.2.3. Phương pháp Picard với sự xấp xỉ liên tục.

Cơ sở của phương pháp Picard là giải chính xác, bởi sự thay thế giá trị y như hàm của x trong phạm vi giá trị x đã cho.

y  g(x)

Đây là biểu thức ước lượng bởi sự thay thế trực tiếp giá trị của x để thu được giá trị tương ứng của y. Cho phương trình vi phân (2.1).

dy = f(x,y)dx

Và tích phân giữa khoảng giới hạn cho x và y.

y 0 y 1 dy = x 0 x 1 f ( x , y ) dx y 0 y 1 dy = x 0 x 1 f ( x , y ) dx size 12{ Int rSub {`y rSub { size 6{0} } } rSup {`y rSub { size 6{1} } } { ital "dy"= Int rSub { size 8{`x rSub { size 6{0} } } } rSup {`x rSub { size 6{1} } } {f \( x,y \) ital "dx"} } } {}

Thì y1y0=x0x1f(x,y)dxy1y0=x0x1f(x,y)dx size 12{y rSub { size 8{1} } - y rSub { size 8{0} } = Int rSub { size 8{`x rSub { size 6{0} } } } rSup {`x rSub { size 6{1} } } {f \( x,y \) ital "dx"} } {}

Hay y1=y0+x0x1f(x,y)dxy1=y0+x0x1f(x,y)dx size 12{y rSub { size 8{1} } =y rSub { size 8{0} } + Int rSub { size 8{`x rSub { size 6{0} } } } rSup {`x rSub { size 6{1} } } {f \( x,y \) ital "dx"} } {} (2.3)

Số hạng tích phân trình bày sự thay đổi trong kết quả của y với sự thay đổi của x từ x0 đến x1. Lời giải có thể thu được bởi sự đánh giá tích phân bằng phương pháp xấp xỉ liên tục.

Ta có thể xem giá trị của y như hàm của x có thể đã thu được bởi sự thay thế y dưới dạng tích phân với y0, cho giá trị ban đầu như sau:

y 1 ( 1 ) = y 0 + x 0 x 1 f ( x , y 0 ) dx y 1 ( 1 ) = y 0 + x 0 x 1 f ( x , y 0 ) dx size 12{y rSub { size 8{1} } rSup { size 8{ \( 1 \) } } =y rSub { size 8{0} } + Int rSub { size 8{`x rSub { size 6{0} } } } rSup {`x rSub { size 6{1} } } {f \( x,y rSub { size 8{0} } \) ital "dx"} } {}

Thực hiện biểu thức tích phân với giá trị mới của y bây giờ được thay thế vào phương trình (2.3) thu được lần xấp xỉ thứ hai cho y như sau:

y 1 ( 2 ) = y 0 + x 0 x 1 f ( x , y 1 ( 1 ) ) dx y 1 ( 2 ) = y 0 + x 0 x 1 f ( x , y 1 ( 1 ) ) dx size 12{y rSub { size 8{1} } rSup { size 8{ \( 2 \) } } =y rSub { size 8{0} } + Int rSub { size 8{`x rSub { size 6{0} } } } rSup {`x rSub { size 6{1} } } {f \( x,y rSub { size 8{1} } rSup { size 8{ \( 1 \) } } \) ` ital "dx"} } {}

Quá trình này có thể lặp lại trong thời gian cần thiết để thu được độ chính xác mong muốn..

Thật vậy, ước lượng tích phân luôn luôn phức tạp thế nhưng phải giả thiết cho biến cố định. Khó khăn và cần thực hiện nhiều lần tích phân, nên đây là mặt hạn chế sự áp dụng của phương pháp này.

Phương pháp Picard có thể áp dụng để giải đồng thời nhiều phương trình như sau:

dy dx = f 1 ( x , y , z ) dy dx = f 1 ( x , y , z ) size 12{ { { ital "dy"} over { ital "dx"} } =f rSub { size 8{1} } \( x,y,z \) } {}

dz dx = f 2 ( x , y , z ) dz dx = f 2 ( x , y , z ) size 12{ { { ital "dz"} over { ital "dx"} } =f rSub { size 8{2} } \( x,y,z \) } {}

Theo công thức, ta có:

y 1 = y 0 + x 0 x 1 f 1 ( x , y 0 , z 0 ) dx y 1 = y 0 + x 0 x 1 f 1 ( x , y 0 , z 0 ) dx size 12{y rSub { size 8{1} } =y rSub { size 8{0} } + Int rSub { size 8{`x rSub { size 6{0} } } } rSup {`x rSub { size 6{1} } } {f rSub { size 8{1} } \( x,y rSub { size 8{0} } ,z rSub { size 8{0} } \) ` ital "dx"} } {}

z 1 = z 0 + x 0 x 1 f 2 ( x , y 0 , z 0 ) dx z 1 = z 0 + x 0 x 1 f 2 ( x , y 0 , z 0 ) dx size 12{z rSub { size 8{1} } =z rSub { size 8{0} } + Int rSub { size 8{`x rSub { size 6{0} } } } rSup {`x rSub { size 6{1} } } {f rSub { size 8{2} } \( x,y rSub { size 8{0} } ,z rSub { size 8{0} } \) ` ital "dx"} } {}

2.2.4. Phương pháp Runge- Kutta.

Trong phương pháp Runge- Kutta sự thay đổi giá trị của biến phụ thuộc là tính toán từ các công thức đã cho, biểu diễn trong điều kiện ước lượng đạo hàm tại những điểm định trước. Từ mỗi giá trị duy nhất chính xác của y cho bởi công thức, phương pháp này không đòi hỏi thay thế lặp lại như phương pháp biến đổi Euler hay tích phân liên tiếp như phương pháp của Picard.

Công thức rút gọn gần đúng xuất phát bởi sự thay thế khai triển chuổi Taylor. Runge- Kutta xấp xỉ bậc hai có thể viết trong công thức.

y1 = y0 + a1k1 + a2k2 (2.4)

Vớik1 = f(x0,y0)h

k2 = f(x0 + b1h, y0 + b2k1)h

Các hệ số a1, a2, b1 và b2 là chính xác. Đầu tiên khai triển f(x0+ b1h, y0+ b2k1) trong chuổi Taylor tại (x0,y0), ta được:

k 2 = f ( x 0 , y 0 ) + b 1 f x 0 h + b 2 k 1 f y 0 + . . . . . h k 2 = f ( x 0 , y 0 ) + b 1 f x 0 h + b 2 k 1 f y 0 + . . . . . h size 12{k rSub { size 8{2} } = left lbrace `f \( x rSub { size 8{0} } ,y rSub { size 8{0} } \) +b rSub { size 8{1} } { { partial f} over { partial x} } \rline rSub { size 8{0} } h+b rSub { size 8{2} } k rSub { size 8{1} } { { partial f} over { partial y} } \rline rSub { size 8{0} } + "." "." "." "." "." right rbrace `h} {}

Thay thế hai điều kiện k1 và k2 vào trong phương trình (2.4), thu được:

y1=y0+(a1+a2)f(x0,y0)h+a2b1fx0h2+a2b2f(x0,y0)fy0h2y1=y0+(a1+a2)f(x0,y0)h+a2b1fx0h2+a2b2f(x0,y0)fy0h2 size 12{y rSub { size 8{1} } =y rSub { size 8{0} } + \( a rSub { size 8{1} } +a rSub { size 8{2} } \) f \( x rSub { size 8{0} } ,y rSub { size 8{0} } \) h+a rSub { size 8{2} } b rSub { size 8{1} } { { partial f} over { partial x} } \rline rSub { size 8{0} } h rSup { size 8{2} } +a rSub { size 8{2} } b rSub { size 8{2} } f \( x rSub { size 8{0} } ,y rSub { size 8{0} } \) { { partial f} over { partial y} } \rline rSub { size 8{0} } h rSup { size 8{2} } } {} (2.5)

Khai triển chuổi Taylor của y tại giá trị (x0,y0) là:

y1=y0+dydx0h+d2ydx20h22+....y1=y0+dydx0h+d2ydx20h22+.... size 12{y rSub { size 8{1} } =y rSub { size 8{0} } + { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } h+ { {d rSup { size 8{2} } y} over { ital "dx" rSup { size 8{2} } } } \rline rSub { size 8{0} } { {h rSup { size 8{2} } } over {2} } + "." "." "." "." } {} (2.6)

Từ dydx0=f(x0,y0)dydx0=f(x0,y0) size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } =f \( x rSub { size 8{0} } ,y rSub { size 8{0} } \) } {}d2ydx20=fx0+fy0f(x0,y0)d2ydx20=fx0+fy0f(x0,y0) size 12{ { {d rSup { size 8{2} } y} over { ital "dx" rSup { size 8{2} } } } \rline rSub { size 8{0} } = { { partial f} over { partial x} } \rline rSub { size 8{0} } + { { partial f} over { partial y} } \rline rSub { size 8{0} } f \( x rSub { size 8{0} } ,y rSub { size 8{0} } \) } {}

Phương trình (2.6) trở thành.

y1=y0+f(x0,y0)h+fx0h22+fy0f(x0,y0)h22......y1=y0+f(x0,y0)h+fx0h22+fy0f(x0,y0)h22...... size 12{y rSub { size 8{1} } =y rSub { size 8{0} } +f \( x rSub { size 8{0} } ,y rSub { size 8{0} } \) h+ { { partial f} over { partial x} } \rline rSub { size 8{0} } { {h rSup { size 8{2} } } over {2} } + { { partial f} over { partial y} } \rline rSub { size 8{0} } f \( x rSub { size 8{0} } ,y rSub { size 8{0} } \) { {h rSup { size 8{2} } } over {2} } "." "." "." "." "." "." } {} (2.7)

Cân bằng các hệ số của phương trình (2.5) và (2.7), ta được:

a1 + a2 =1; a2b1 = 1/2; a2b2 = 1/2.

Chọn giá trị tùy ý cho a1

a1 = 1/2

Thì a2 = 1/2; b1 = 1; b2 = 1.

Thay thế giá trị này vào trong phương trình (2.4), công thức gần đúng bậc hai Runge-Kutta là:

y 1 = y 0 + k 1 + k 2 y 1 = y 0 + k 1 + k 2 size 12{y rSub { size 8{1} } =y rSub { size 8{0} } + {1} wideslash {2} k rSub { size 8{1} } + {1} wideslash {2} k rSub { size 8{2} } } {}

Vớik1 = f(x0,y0)h

k2 = f(x0+ h, y0 + k1)h

Vì thế.

Δy = ( k 1 + k 2 ) Δy = ( k 1 + k 2 ) size 12{Δy= {1} wideslash {2} \( k rSub { size 8{1} } +k rSub { size 8{2} } \) } {}

Áp dụng của phương pháp Runge-Kutta cho việc xấp xỉ bậc hai đòi hỏi sự tính toán của k1 và k2­. Sai số trong lần xấp xỉ là bậc h3 bởi vì chuổi đã cắt sau điều kiện bậc hai.

Tông quát công thức xấp xỉ bậc bốn Runge-Kutta là:

y1=y0+a1k1+a2k2+a3k3+a4k4y1=y0+a1k1+a2k2+a3k3+a4k4 size 12{y rSub { size 8{1} } =y rSub { size 8{0} } +a rSub { size 8{1} } k rSub { size 8{1} } +a rSub { size 8{2} } k rSub { size 8{2} } +a rSub { size 8{3} } k rSub { size 8{3} } +a rSub { size 8{4} } k rSub { size 8{4} } } {} (2.8)

Với k1 = f(x0,y0)h

k2 = f(x0 + b1h, y0 + b2k1)h

k3 = f(x0 + b3h, y0 + b4k2)h

k4 = f(x0 + b5h, y0 + b6k3)h

Tiếp theo thủ tục giống như dùng cho lần xấp xỉ bậc hai, hệ số trong phương trình (2.8) thu được là:

a1 = 1/6; a2 = 2/6; a3 = 2/6; a4 = 1/6.

Và b1 = 1/2; b2 = 1/2; b3 = 1/2; b4 = 1/2; b5 = 1; b6 = 1.

Thay thế các giá trị vào trong phương trình (2.8), phương trình xấp xỉ bậc bốn Runge-Kutta trở thành.

y 1 = y 0 + ( k 1 + 2k 2 + 2k 3 + k 4 ) y 1 = y 0 + ( k 1 + 2k 2 + 2k 3 + k 4 ) size 12{y rSub { size 8{1} } =y rSub { size 8{0} } + {1} wideslash {6} \( k rSub { size 8{1} } +2k rSub { size 8{2} } +2k rSub { size 8{3} } +k rSub { size 8{4} } \) } {}

Vớik1 = f(x0,y0)h

k 2 = f ( x 0 + h 2 , y 0 + k 1 2 ) h k 2 = f ( x 0 + h 2 , y 0 + k 1 2 ) h size 12{k rSub { size 8{2} } =f \( x rSub { size 8{0} } + { {h} over {2} } ,y rSub { size 8{0} } + { {k rSub { size 8{1} } } over {2} } \) h} {}

k 3 = f ( x 0 + h 2 , y 0 + k 2 2 ) h k 3 = f ( x 0 + h 2 , y 0 + k 2 2 ) h size 12{k rSub { size 8{3} } =f \( x rSub { size 8{0} } + { {h} over {2} } ,y rSub { size 8{0} } + { {k rSub { size 8{2} } } over {2} } \) h} {}

k 4 = f ( x 0 + h , y 0 + k 3 ) h k 4 = f ( x 0 + h , y 0 + k 3 ) h size 12{k rSub { size 8{4} } =f \( x rSub { size 8{0} } +h,y rSub { size 8{0} } +k rSub { size 8{3} } \) h} {}

Như vậy, sự tính toán của y theo công thức đòi hỏi sự tính toán các giá trị của k1, k2, k3 và k4 :

y = 1/6(k1+2k2+2k3+k4)

Sai số trong sự xấp xỉ là bậc h5.

Công thức xấp xỉ bậc bốn Runge-Kutta cho phép giải đồng thời nhiều phương trình vi phân.

dy dx = f ( x , y , z ) dy dx = f ( x , y , z ) size 12{ { { ital "dy"} over { ital "dx"} } =f \( x,y,z \) } {}

dz dx = g ( x , y , z ) dz dx = g ( x , y , z ) size 12{ { { ital "dz"} over { ital "dx"} } =g \( x,y,z \) } {}

Ta co:

y­1­ = y0+1/6 (k1+2k2+2k3+k4)

z­1­ = z0+1/6 (l1+2l2+2l3+l4)

Với:k1= f(x0,y0,z0)h

k 2 = f ( x 0 + h 2 , y 0 + k 1 2 z 0 + l 1 2 ) h k 2 = f ( x 0 + h 2 , y 0 + k 1 2 z 0 + l 1 2 ) h size 12{k rSub { size 8{2} } =f \( x rSub { size 8{0} } + { {h} over {2} } ,y rSub { size 8{0} } + { {k rSub { size 8{1} } } over {2} } z rSub { size 8{0} } + { {l rSub { size 8{1} } } over {2} } \) h} {}

k 3 = f ( x 0 + h 2 , y 0 + k 2 2 z 0 + l 2 2 ) h k 3 = f ( x 0 + h 2 , y 0 + k 2 2 z 0 + l 2 2 ) h size 12{k rSub { size 8{3} } =f \( x rSub { size 8{0} } + { {h} over {2} } ,y rSub { size 8{0} } + { {k rSub { size 8{2} } } over {2} } z rSub { size 8{0} } + { {l rSub { size 8{2} } } over {2} } \) h} {}

k4 = f(x0 + h, y0 + k3,z0 + l3)h

l1 = g(x0,y0,z0)h

l 2 = g ( x 0 + h 2 , y 0 + k 1 2 z 0 + l 1 2 ) h l 2 = g ( x 0 + h 2 , y 0 + k 1 2 z 0 + l 1 2 ) h size 12{l rSub { size 8{2} } =g \( x rSub { size 8{0} } + { {h} over {2} } ,y rSub { size 8{0} } + { {k rSub { size 8{1} } } over {2} } z rSub { size 8{0} } + { {l rSub { size 8{1} } } over {2} } \) h} {}

l 3 = g ( x 0 + h 2 , y 0 + k 2 2 z 0 + l 2 2 ) h l 3 = g ( x 0 + h 2 , y 0 + k 2 2 z 0 + l 2 2 ) h size 12{l rSub { size 8{3} } =g \( x rSub { size 8{0} } + { {h} over {2} } ,y rSub { size 8{0} } + { {k rSub { size 8{2} } } over {2} } z rSub { size 8{0} } + { {l rSub { size 8{2} } } over {2} } \) h} {}

l4 = g(x0 + h, y0 + k3,z0 + l3)h

2.2.5. Phương pháp dự đoán sửa đổi.

Phương pháp dựa trên cơ sở ngoại suy, hay tích phân vượt trước, và lặp lại nhiều lần việc giải phương trình vi phân.

dydx=f(x,y)dydx=f(x,y) size 12{ { { ital "dy"} over { ital "dx"} } =f \( x,y \) } {} (2.9)

Được gọi là phương pháp dự đoán sửa đổi. Thủ tục cơ bản trong phương pháp dự đoán sửa đổi là xuất phát từ điểm (xn,yn) đến điểm (xn+1, yn+1). Thì thu được dydxn+1dydxn+1 size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{n+1} } } {}từ phương trình vi phân và sửa đổi giá trị yn+1 xấp xỉ công thức chính xác.

Loại đơn giản của công thức dự đoán phương pháp của Euler là:

yn+1 = yn + yn’h (2.10)

Với: yn'=dydxnyn'=dydxn size 12{y rSub { size 8{n} } rSup { size 8{'} } = { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{n} } } {}

Công thức chính xác không dùng trong phương pháp Euler. Mặc dù, trong phương pháp biến đổi Euler giá trị gần đúng của yn+1 thu được từ công thức dự đoán (2.10) và giá trị thay thế trong phương trình vi phân (2.9) chính là y’n+1. Thì giá trị chính xác cho yn+1 thu được từ công thức biến đổi của phương pháp là:

yn+1=yn+(y'n+1+y'n)h2yn+1=yn+(y'n+1+y'n)h2 size 12{y rSub { size 8{n+1} } =y rSub { size 8{n} } + \( y' rSub { size 8{n+1} } +y' rSub { size 8{n} } \) { {h} over {2} } } {} (2.11)

Giá trị thay thế trong phương trình vi phân (2.9) thu được có sự đánh giá chính xác hơn cho y’n+1, nó luôn luôn thay thế trong phương trình (2.11) làm cho yn+1 chính xác hơn. Quá trình tiếp tục lặp lại cho đến khi hai giá trị tính toán liên tiếp của yn+1 từ phương trình (2.11) trùng với giá trị mong muốn chấp nhận được.

Phương pháp dự đoán biến đổi kinh điển của Milne. Dự đoán của Milne và công thức biến đổi, theo ông là:

y n + 1 ( 0 ) = y n 3 + 4h 3 ( 2y ' n 2 y ' n 1 + 2y ' n ) y n + 1 ( 0 ) = y n 3 + 4h 3 ( 2y ' n 2 y ' n 1 + 2y ' n ) size 12{y rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } =y rSub { size 8{n - 3} } + { {4h} over {3} } \( 2y' rSub { size 8{n - 2} } - y' rSub { size 8{n - 1} } +2y' rSub { size 8{n} } \) } {}

yn+1=yn1+h3(y'n1+4y'n+y'n+1)yn+1=yn1+h3(y'n1+4y'n+y'n+1) size 12{y rSub { size 8{n+1} } =y rSub { size 8{n - 1} } + { {h} over {3} } \( y' rSub { size 8{n - 1} } +4y' rSub { size 8{n} } +y' rSub { size 8{n+1} } \) } {}

Với: y'n+1=f(xn+1,yn+1(0))y'n+1=f(xn+1,yn+1(0)) size 12{y' rSub { size 8{n+1} } =f \( x rSub { size 8{n+1} } ,y rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } \) } {}

Bắt đầu của sự tính toán đòi hỏi biết bốn giá trị của y. Có thể đã tính toán bởi Runge-Kutta hay một số phương pháp số trước khi sử dụng công thức dự đoán sửa đổi của Milne. Sai số trong phương pháp là bậc h5.

Trong trường hợp tổng quát, phương pháp mong muốn chọn h đủ nhỏ nên chỉ vài lần lặp là đòi hỏi thu được yn+1 hoàn toàn chính xác như mong muốn.

Phương pháp có thể mở rộng cho phép giải một số phương trình vi phân đồng thời. Phương pháp dự đoán sửa đổi là áp dụng độc lập đối với mỗi phương trình vi phân như một phương trình vi phân đơn giản. Vì vậy, thay thế giá trị cho tất cả các biến phụ thuộc vào trong mỗi phương trình vi phân là đòi hỏi sự đánh giá đạo hàm tại (xn+1, yn+1).

2.3. GIẢI PHƯƠNG TRÌNH VI PHÂN BẬC CAO.

Trong kỹ thuật trước đây mô tả cho việc giải phương trình vi phân bậc nhất cũng có thể áp dụng cho việc giải phương trình vi phân bậc cao bằng sự đưa vào của biến phụ. Ví dụ, cho phương trình vi phân bậc hai.

a d 2 y dx 2 + b dy dx + cy = 0 a d 2 y dx 2 + b dy dx + cy = 0 size 12{a { {d rSup { size 8{2} } y} over { ital "dx" rSup { size 8{2} } } } +b { { ital "dy"} over { ital "dx"} } + ital "cy"=0} {}

Với điều kiện ban đầu x0, y0, và dydx0dydx0 size 12{ { { ital "dy"} over { ital "dx"} } \rline rSub { size 8{0} } } {}thì phương trình có thể được viết lại như hai phương trình vi phân bậc nhất.

dy dx = y ' dy dx = y ' size 12{ { { ital "dy"} over { ital "dx"} } =y'} {}

d 2 y dx 2 = dy ' dx = by ' + cy a d 2 y dx 2 = dy ' dx = by ' + cy a size 12{ { {d rSup { size 8{2} } y} over { ital "dx" rSup { size 8{2} } } } = { { ital "dy"'} over { ital "dx"} } = - { { ital "by"'+ ital "cy"} over {a} } } {}

Một trong những phương pháp mô tả trước đây có thể là việc làm đi tìm lời giải cho hai phương trình vi phân bậc nhất đồng thời.

Theo cách tương tự, một vài phương trình hay hệ phương trình bậc cao có thể quy về hệ phương trình vi phân bậc nhất.

2.4. VÍ DỤ VỀ GIẢI PHƯƠNG TRÌNH VI PHÂN BẰNG PHƯƠNG PHÁP SỐ.

Giải phương trình vi phân sẽ minh họa bằng sự tính toán dòng điện cho mạch RL nối tiếp.

e(t)LRt = 0i(t)Hình 2.4: Sự biểu diễn của mạch điện RLCho mạch điện RL trong hình 2.4 sức điện động hiệu dụng khi đóng khóa là:

e(t) = 5t 0  t  0,2

e(t) = 1 t > 0,2

Điện trở cho theo đơn vị ohms là.

R = 1+3i2

Và điện cảm theo đơn vị henrys là.

L = 1

Tìm dòng điện trong mạch điện theo các phương pháp sau:

  1. Euler’s
  2. Biến đổi Euler.
  3. Xấp xỉ bậc bốn Runge-Kutta
  4. Milne’s
  5. Picard’s

Bài giải:

Phương trình vi phân của mạch điện là.

L di dt + Ri = e ( t ) L di dt + Ri = e ( t ) size 12{L { { ital "di"} over { ital "dt"} } + ital "Ri"=e \( t \) } {}

Thay thế cho R và L ta có:

di dt + ( 1 + 3i 2 ) i = e ( t ) di dt + ( 1 + 3i 2 ) i = e ( t ) size 12{ { { ital "di"} over { ital "dt"} } + \( 1+3i rSup { size 8{2} } \) i=e \( t \) } {}

Điều kiện ban đầu tại t = 0 thì e0 = 0 và i0 = 0. Khoảng chọn cho biến độc lập là:

t = 0,025.

a. Phương trình theo phương pháp Euler là.

Δi n = di dt n Δt Δi n = di dt n Δt size 12{Δi rSub { size 8{n} } = { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n} } Δt} {}

in+1 = in +in

Với didtn=en(1+3in2)indidtn=en(1+3in2)in size 12{ { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n} } =e rSub { size 8{n} } - \( 1+3i rSub { size 8{n} } rSup { size 8{2} } \) i rSub { size 8{n} } } {}

Thay thế giá trị ban đầu vào trong phương trình vi phân, dydt0=0dydt0=0 size 12{ { { ital "dy"} over { ital "dt"} } \rline rSub { size 8{0} } =0} {} và i0. Vì thế, dòng điện i1 = 0. Tại t1 = 0,025; e1 = 0,125 và didt1=0,125{1+3(0)2}0=0,125didt1=0,125{1+3(0)2}0=0,125 size 12{ { { ital "di"} over { ital "dt"} } \rline rSub { size 8{1} } =0,"125" - lbrace 1+3 \( 0 \) rSup { size 8{2} } rbrace 0=0,"125"} {}

i1 = (0,125)0,025 = 0,00313

Thì

i2 = 0 + 0,00313 = 0,00313

Lập bảng kê kết quả lời giải đưa vào trong bảng 2.1

Bảng 2.1: Giải bằng phương pháp Euler

Bảng 1
in=in1+didtn1Δtin=in1+didtn1Δt size 12{i rSub { size 8{n} } =i rSub { size 8{n - 1} } + { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n - 1} } Δt} {}didtn=en(1+3in2)indidtn=en(1+3in2)in size 12{ { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n} } =e rSub { size 8{n} } - \( 1+3i rSub { size 8{n} } rSup { size 8{2} } \) i rSub { size 8{n} } } {}n Thời giantn Sức điện độngen Dòng  
0123456789101112 0,0000,0250,0500,0750,1000,1250,1500,1750,2000,2250,2500,2750,300 0,0000,1250,2500,2500,3750,5000.6250,7500,8751,0001,0001,0001,000 0,000000,000000,003130,009300,018440,030480,45340,062950,083230,106110,128370,150000,17100 0,000000,125000,246870,365700,481540,594440,704380,811300,915040,890310,865280,83988

b. Phương trình của phương pháp biến đổi Euler là.

Δi n ( 0 ) = di dt n Δt Δi n ( 0 ) = di dt n Δt size 12{Δi rSub { size 8{n} } rSup { size 8{ \( 0 \) } } = { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n} } Δt} {}

i n + 1 ( 0 ) = i n + Δi n ( 0 ) i n + 1 ( 0 ) = i n + Δi n ( 0 ) size 12{i rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } =i rSub { size 8{n} } +Δi rSub { size 8{n} } rSup { size 8{ \( 0 \) } } } {}

Δi n ( 1 ) = di dt n + di dt n + 1 ( 0 ) 2 Δt Δi n ( 1 ) = di dt n + di dt n + 1 ( 0 ) 2 Δt size 12{Δi rSub { size 8{n} } rSup { size 8{ \( 1 \) } } = left ( { { { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n} } + { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } } over {2} } right )Δt} {}

i n + 1 ( 1 ) = i n + Δi n ( 1 ) i n + 1 ( 1 ) = i n + Δi n ( 1 ) size 12{i rSub { size 8{n+1} } rSup { size 8{ \( 1 \) } } =i rSub { size 8{n} } +Δi rSub { size 8{n} } rSup { size 8{ \( 1 \) } } } {}

Với didtn+1(0)=en+1{1+3(in+1(0))2}in+1(0)didtn+1(0)=en+1{1+3(in+1(0))2}in+1(0) size 12{ { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } =e rSub { size 8{n+1} } - lbrace 1+3 \( i rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } \) rSup { size 8{2} } rbrace i rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } } {}

Thay thế giá trị ban đầu e0 = 0 và i0 = 0 vào trong phương trình vi phân didx0=0didx0=0 size 12{ { { ital "di"} over { ital "dx"} } \rline rSub { size 8{0} } =0} {}

Do đó: Δi0(0)=0Δi0(0)=0 size 12{Δi rSub { size 8{0} } rSup { size 8{ \( 0 \) } } =0} {}; i1(0)=0i1(0)=0 size 12{i rSub { size 8{1} } rSup { size 8{ \( 0 \) } } =0} {}.

Thay thế vào trong phương trình vi phân i1(0)=0i1(0)=0 size 12{i rSub { size 8{1} } rSup { size 8{ \( 0 \) } } =0} {} và e1 = 0,125

di dt 1 ( 0 ) = 0, 125 { 1 + 3 ( 0 ) 2 } 0 = 0, 125 di dt 1 ( 0 ) = 0, 125 { 1 + 3 ( 0 ) 2 } 0 = 0, 125 size 12{ { { ital "di"} over { ital "dt"} } \rline rSub { size 8{1} } rSup { size 8{ \( 0 \) } } =0,"125" - lbrace 1+3 \( 0 \) rSup { size 8{2} } rbrace 0=0,"125"} {}

Δi0(1)=(0,125+02)0,025=0,00156Δi0(1)=(0,125+02)0,025=0,00156 size 12{Δi rSub { size 8{0} } rSup { size 8{ \( 1 \) } } = \( { {0,"125"+0} over {2} } \) 0,"025"=0,"00156"} {}

Nên

i 1 ( 1 ) = 0 + 0, 00156 = 0, 00156 i 1 ( 1 ) = 0 + 0, 00156 = 0, 00156 size 12{i rSub { size 8{1} } rSup { size 8{ \( 1 \) } } =0+0,"00156"=0,"00156"} {}

Trong lời giải ví dụ cho phương pháp, không thực hiện lặp lại in+1(1)=in+1in+1(1)=in+1 size 12{i rSub { size 8{n+1} } rSup { size 8{ \( 1 \) } } =i rSub { size 8{n+1} } } {}. Bài giải thu được bằng phương pháp biến đổi Euler được đưa vào trong bảng 2.2.

didtndidtn size 12{ { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n} } } {}didtn+1(0)didtn+1(0) size 12{ { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } } {}Bảng 2.2: Bài giải bằng phương pháp biến đổi Euler.

Bảng 2
Δin(0)Δin(0) size 12{Δi rSub { size 8{n} } rSup { size 8{ \( 0 \) } } } {}en+1en+1 size 12{e rSub { size 8{n+1} } } {}in+1(0)in+1(0) size 12{i rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } } {}Δin(1)Δin(1) size 12{Δi rSub { size 8{n} } rSup { size 8{ \( 1 \) } } } {}n Thời Sức Dòng Gian điện điện intn động en
0123456789101112 0,000 0,000 0,00000 0,00000 0,00000 0,125 0,00000 0,12500 0,001560,025 0,125 0,00156 0,12344 0,00309 0,250 0,00465 0,24535 0,004610,050 0,250 0,00617 0,34383 0,00610 0,375 0,01227 0,36272 0,007580,075 0,375 0,01375 0,36124 0,00903 0,500 0,02278 0,47718 0,01048 FIXME: A LIST CAN NOT BE A TABLE ENTRY. 0,500 0,02423 0,47573 0,01189 0,625 0,03612 0,58874 0,01331 FIXME: A LIST CAN NOT BE A TABLE ENTRY. 0,625 0,03754 0,58730 0,01468 0,750 0,05222 0,69735 0,01606 FIXME: A LIST CAN NOT BE A TABLE ENTRY. 0,750 0,05360 0,69594 0,01740 0,875 0,07100 0,80293 0,018740,175 0,875 0,07234 0,80152 0,02004 1,000 0,09238 0,90525 0,021330,200 1,000 0,09367 0,90386 0,02260 1,000 0,11627 0,87901 0,022290,225 1,000 0,11596 0,87936 0,02198 1,000 0,13794 0,85419 0,021670,250 1,000 0,13763 0,85455 0,02136 1,000 0,15899 0,82895 0,021040,275 1,000 0,15867 0,82935 0,02073 1,000 0,17940 0,80328 0,020410,300 1,000 0,17908

c. Phương trình dùng phương pháp Runge-Kutta để giải.

di dt = e ( t ) ( 1 + 3i 2 ) i di dt = e ( t ) ( 1 + 3i 2 ) i size 12{ { { ital "di"} over { ital "dt"} } =e \( t \) - \( 1+3i rSup { size 8{2} } \) i} {}

Ta có:

k 1 = { e ( t n ) ( 1 + 3i n 2 ) i n } Δt k 1 = { e ( t n ) ( 1 + 3i n 2 ) i n } Δt size 12{k rSub { size 8{1} } = lbrace e \( t rSub { size 8{n} } \) - \( 1+3i rSub { size 8{n} } rSup { size 8{2} } \) i rSub { size 8{n} } rbrace Δt} {}

k 2 = e ( t n + Δt 2 ) 1 + 3 i n + k 1 2 2 . i n + k 1 2 Δt k 2 = e ( t n + Δt 2 ) 1 + 3 i n + k 1 2 2 . i n + k 1 2 Δt size 12{k rSub { size 8{2} } = left lbrace e \( t rSub { size 8{n} } + { {Δt} over {2} } \) - left [1+3 left (i rSub { size 8{n} } + { {k rSub { size 8{1} } } over {2} } right ) rSup { size 8{2} } right ]` "." ` left (i rSub { size 8{n} } + { {k rSub { size 8{1} } } over {2} } right ) right rbrace Δt} {}

k 3 = e ( t n + Δt 2 ) 1 + 3 i n + k 2 2 2 . i n + k 2 2 Δt k 3 = e ( t n + Δt 2 ) 1 + 3 i n + k 2 2 2 . i n + k 2 2 Δt size 12{k rSub { size 8{3} } = left lbrace e \( t rSub { size 8{n} } + { {Δt} over {2} } \) - left [1+3 left (i rSub { size 8{n} } + { {k rSub { size 8{2} } } over {2} } right ) rSup { size 8{2} } right ]` "." ` left (i rSub { size 8{n} } + { {k rSub { size 8{2} } } over {2} } right ) right rbrace Δt} {}

k 4 = { e ( t n + Δt ) 1 + 3 ( i n + k 3 ) 2 . ( i n + k 3 ) } Δt k 4 = { e ( t n + Δt ) 1 + 3 ( i n + k 3 ) 2 . ( i n + k 3 ) } Δt size 12{k rSub { size 8{4} } = lbrace e \( t rSub { size 8{n} } +Δt \) - left [1+3 \( i rSub { size 8{n} } +k rSub { size 8{3} } \) rSup { size 8{2} } right ]` "." ` \( i rSub { size 8{n} } +k rSub { size 8{3} } \) rbrace Δt} {}

Δi n = ( k 1 + 2k 2 + 2k 3 + k 4 ) Δi n = ( k 1 + 2k 2 + 2k 3 + k 4 ) size 12{Δi rSub { size 8{n} } = {1} wideslash {6} \( k rSub { size 8{1} } +2k rSub { size 8{2} } +2k rSub { size 8{3} } +k rSub { size 8{4} } \) } {}

in+1 = in + in

Với:

e(tn) = en

e ( t n + Δt 2 ) = e n + e n + 1 2 e ( t n + Δt 2 ) = e n + e n + 1 2 size 12{e \( t rSub { size 8{n} } + { {Δt} over {2} } \) = { {e rSub { size 8{n} } +e rSub { size 8{n+1} } } over {2} } } {}

e(tn + t) = en+1

Thay thế giá trị ban đầu tìm được k1:

k1 = 0.

Tìm được k2­:

k 2 = 0 + 0, 125 2 1 + 3 ( 0 ) 2 0 0, 025 = 0, 00156 k 2 = 0 + 0, 125 2 1 + 3 ( 0 ) 2 0 0, 025 = 0, 00156 size 12{k rSub { size 8{2} } = left lbrace { {0+0,"125"} over {2} } - left [1+3 \( 0 \) rSup { size 8{2} } right ]`0 right rbrace 0,"025"=0,"00156"} {}

Tìm được k3:

k 3 = 0 + 0, 125 2 1 + 3 0, 00156 2 2 0, 00156 2 0, 025 = 0, 00154 k 3 = 0 + 0, 125 2 1 + 3 0, 00156 2 2 0, 00156 2 0, 025 = 0, 00154 size 12{k rSub { size 8{3} } = left lbrace { {0+0,"125"} over {2} } - left [1+3 left ( { {0,"00156"} over {2} } right ) rSup { size 8{2} } right ] { {0,"00156"} over {2} } right rbrace 0,"025"=0,"00154"} {}

Tìm được k4:

­ k4=0+0,1251+3(0,00154)20,001540,025=0,00309k4=0+0,1251+3(0,00154)20,001540,025=0,00309 size 12{k rSub { size 8{4} } = left lbrace 0+0,"125" - left [1+3 \( 0,"00154" \) rSup { size 8{2} } right ]`0,"00154" right rbrace `0,"025"=0,"00309"} {}

Thì

Δi 0 = ( 0 + 0, 00312 + 0, 00308 + 0, 00309 ) = 0, 00155 Δi 0 = ( 0 + 0, 00312 + 0, 00308 + 0, 00309 ) = 0, 00155 size 12{Δi rSub { size 8{0} } = {1} wideslash {6} \( 0+0,"00312"+0,"00308"+0,"00309" \) =0,"00155"} {}

Và i1 = i0 + i0 = 0+ 0,00155 = 0,00155

Bài giải thu được bằng phương pháp Runge-Kutta được đưa vào trong bảng 2.3.

d. Công thức dự đoán sửa đổi của phương pháp Milne là.

i n + 1 ( 0 ) = i n 3 + 4Δt 3 ( 2i ' n 2 i ' n 1 + 2i ' n ) i n + 1 ( 0 ) = i n 3 + 4Δt 3 ( 2i ' n 2 i ' n 1 + 2i ' n ) size 12{i rSub { size 8{n+1} } rSup { size 8{ \( 0 \) } } =i rSub { size 8{n - 3} } + { {4Δt} over {3} } \( 2i' rSub { size 8{n - 2} } - i' rSub { size 8{n - 1} } +2i' rSub { size 8{n} } \) } {}

i n + 1 = i n 1 + Δt 3 ( i ' n 1 + 4i ' n + i ' n + 1 ) i n + 1 = i n 1 + Δt 3 ( i ' n 1 + 4i ' n + i ' n + 1 ) size 12{i rSub { size 8{n+1} } =i rSub { size 8{n - 1} } + { {Δt} over {3} } \( i' rSub { size 8{n - 1} } +4i' rSub { size 8{n} } +i' rSub { size 8{n+1} } \) } {}

Với

i ' n = di dt n i ' n = di dt n size 12{i' rSub { size 8{n} } = { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n} } } {}

di dt n = e n ( 1 + 3i n 2 ) i n di dt n = e n ( 1 + 3i n 2 ) i n size 12{ { { ital "di"} over { ital "dt"} } \rline rSub { size 8{n} } =e rSub { size 8{n} } - \( 1+3i rSub { size 8{n} } rSup { size 8{2} } \) i rSub { size 8{n} } } {}

Các giá trị ban đầu đòi hỏi phải thu được từ lời giải của phương pháp Runge-Kutta.

Với i0 = 0; i1 = 0,00155; i2 = 0,00615; i3 = 0,01372.

Thay thế vào phương trình vi phân, ta có:

i’0 = 0; i’1 = 0,12345; i’2 = 0,23485; i’3 = 0,36127.

Bắt đầu tại t4 = 0,100 và thay thế vào trong công thức dự đoán, ước lượng đầu tiên cho i4 là:

i 4 ( 0 ) = 0 + ( 0, 025 ) 2 ( 0, 12345 ) 0, 24385 + 2 ( 0, 36127 ) = 0, 02418 i 4 ( 0 ) = 0 + ( 0, 025 ) 2 ( 0, 12345 ) 0, 24385 + 2 ( 0, 36127 ) = 0, 02418 size 12{i rSub { size 8{4} } rSup { size 8{ \( 0 \) } } =0+ {4} wideslash {3} \( 0,"025" \) left [2 \( 0,"12345" \) - 0,"24385"+2 \( 0,"36127" \) right ]=0,"02418"} {}

Thay thế e4 = 0,500 và i4 = 0,02418 vào trong phương trình vi phân, ta được:

i’4 = 0,500 [ 1 + 3(0,02418)2]0,02418 = 0,47578

Dự đoán và giá trị chính xác, chỉ khác nhau một số hàng thập phân vì vậy không đòi hỏi lặp lại nhiều lần. Kết quả sau từng bước được ghi vào bảng 2.4. Tại t9 giá trị dự đoán của dòng điện là 0,11742 nhưng trong khi giá trị chính xác là 0,11639. Việc thực hiện lặp lại bởi sự thay thế giá trị chính xác trong phương trình vi phân đã thu được i’9 = 0,87888. Cứ lần lượt dùng trong công thức sửa đổi để thu được ước lượng thứ hai cho i9 = 0,11640, trước khi kiểm tra giá trị chính xác. Thực hiện lặp lại trong tất cả các bước để đảm bảo yêu cầu chính xác.

Bảng 3
Bảng 2.3: Giải bằng phương pháp Runge-Kutta
Thời Sức Dòng en+ en+1 k1 k2 gian điện điện k1 -------- in + --- k2 in + --- k3 en+1 in + k3 k4 intn động in 2 2 2 en 0,000 0,000 0,00000 0,00000 0,0625 0,00000 0,00156 0,00078 0,00154 0,125 0,00154 0,00309 0,001550,025 0,125 0,00155 0,00309 0,1875 0,00310 0,00461 0,00386 0,00459 0,250 0,00614 0,00610 0,004600,050 0,250 0,00615 0,00610 0,3125 0,00920 0,00758 0,00994 0,00756 0,375 0,01371 0,00903 0,007570,075 0,375 0,01372 0,00903 0,4375 0,01824 0,01048 0,01896 0,01046 0,500 0,02418 0,01189 0,010470,100 0,500 0,02419 0,01189 0,5625 0,03014 0,01331 0,03084 0,01329 0,625 0,03748 0,01468 0,013300,125 0,625 0,03749 0,01468 0,6875 0,04483 0,01606 0,04552 0,01604 0,750 0,05353 0,01740 0,01605 FIXME: A LIST CAN NOT BE A TABLE ENTRY. 0,750 0,05354 0,01740 0,8125 0,06224 0,01874 0,06291 0,01872 0,875 0,07226 0,02004 0,018730,175 0,875 0,07227 0,02004 0,9375 0,08229 0,02134 0,08294 0,02132 1,000 0,09359 0,02260 0,021330,200 1,000 0,09360 0,02260 1,0000 0,10490 0,02229 0,10475 0,02230 1,000 0,11590 0,02199 0,022300,225 1,000 0,11590 0,02199 1,0000 0,12690 0,02167 0,12674 0,02168 1,000 0,13758 0,02137 0,021680,250 1,000 0,13758 0,02137 1,0000 0,14827 0,02105 0,14811 0,02105 1,000 0,15863 0,02073 0,021050,275 1,000 0,15863 0,02073 1,0000 0,16900 0,02041 0,16884 0,02042 1,000 0,17905 0,02009 0,02041
n 0123456789101112

Bảng 2.4: Bài giải bằng phương pháp của Milne.

Bảng 4
N Thời gian Sức điện Dòng điện Dòng điệntn động en (dự đoán) in i’n (sửa đổi) in
456789101112 0,100 0,500 0,02418 0,47578 0,024190,125 0,625 0,03748 0,58736 0,037480,150 0,750 0,05353 0,69601 0,053530,175 0,875 0,07226 0,80161 0,072260,200 1,000 0,09359 0,90395 0,093580,225 1,000 0,11742 0,87772 0,116390,87888 0,11640+0,250 1,000 0,13543 0,85712 0,137550,85464 0,13753+0,275 1,000 0,16021 0,82745 0,159110,82881 0,15912+0,300 1,000 0,17894 0,80387 0,178980,80382 0,17898+
 

+ : giá trị sửa đổi thứ hai thu được bởi vòng lặp

d. Phương trình dùng phương pháp Picard hàm tương đương khởi đầu cho i, cận i0 = 0 là:

i = i 0 + 0 t e ( t ) i 3i 3 dt i = i 0 + 0 t e ( t ) i 3i 3 dt size 12{i=i rSub { size 8{0} } + Int rSub { size 8{`0} } rSup { size 8{`t} } { left [e \( t \) - i - 3i rSup { size 8{3} } right ]} ` ital "dt"} {}

Thay thế e(t) = 5t và giá trị ban đầu i0 = 0

i ( 1 ) = 0 t 5 t dt = 5t 2 2 i ( 1 ) = 0 t 5 t dt = 5t 2 2 size 12{i rSup { size 8{ \( 1 \) } } = Int rSub { size 8{`0} } rSup { size 8{`t} } {5`t` ital "dt"= { {5t rSup { size 8{2} } } over {2} } } } {}

Thay i(1) cho i trong phương trình tích phân, thu được:

i ( 2 ) = 0 t 5t 5t 2 2 375 t 6 8 dt = 5t 2 2 5t 3 6 375 t 7 56 i ( 2 ) = 0 t 5t 5t 2 2 375 t 6 8 dt = 5t 2 2 5t 3 6 375 t 7 56 size 12{i rSup { size 8{ \( 2 \) } } = Int rSub { size 8{`0} } rSup { size 8{`t} } { left (5t - { {5t rSup { size 8{2} } } over {2} } - { {"375"t rSup { size 8{6} } } over {8} } right )} ` ital "dt"= { {5t rSup { size 8{2} } } over {2} } - { {5t rSup { size 8{3} } } over {6} } - { {"375"t rSup { size 8{7} } } over {"56"} } } {}

Quá trình tiếp tục, ta được:

i ( 3 ) = 0 t 5t 5t 2 2 + 5t 3 6 375 t 6 8 + 375 t 7 7 125 t 8 8 + . . . . dt i ( 3 ) = 0 t 5t 5t 2 2 + 5t 3 6 375 t 6 8 + 375 t 7 7 125 t 8 8 + . . . . dt size 12{i rSup { size 8{ \( 3 \) } } = Int rSub { size 8{`0} } rSup { size 8{`t} } { left (5t - { {5t rSup { size 8{2} } } over {2} } + { {5t rSup { size 8{3} } } over {6} } - { {"375"t rSup { size 8{6} } } over {8} } + { {"375"t rSup { size 8{7} } } over {7} } - { {"125"t rSup { size 8{8} } } over {8} } + "." "." "." "." right )} ` ital "dt"} {}

= 5t 2 2 5t 3 6 + 5t 4 24 375 t 7 56 + . . . . = 5t 2 2 5t 3 6 + 5t 4 24 375 t 7 56 + . . . . size 12{ {}= { {5t rSup { size 8{2} } } over {2} } - { {5t rSup { size 8{3} } } over {6} } + { {5t rSup { size 8{4} } } over {"24"} } - { {"375"t rSup { size 8{7} } } over {"56"} } + "." "." "." "." } {}

i ( 4 ) = 0 t 5t 5t 2 2 + 5t 3 6 5t 4 24 375 t 6 8 + 375 t 7 7 + . . . . dt i ( 4 ) = 0 t 5t 5t 2 2 + 5t 3 6 5t 4 24 375 t 6 8 + 375 t 7 7 + . . . . dt size 12{i rSup { size 8{ \( 4 \) } } = Int rSub { size 8{`0} } rSup { size 8{`t} } { left (5t - { {5t rSup { size 8{2} } } over {2} } + { {5t rSup { size 8{3} } } over {6} } - { {5t rSup { size 8{4} } } over {"24"} } - { {"375"t rSup { size 8{6} } } over {8} } + { {"375"t rSup { size 8{7} } } over {7} } + "." "." "." "." right )} ` ital "dt"} {}

= 5t 2 2 5t 3 6 + 5t 4 24 t 5 24 375 t 7 56 + . . . . = 5t 2 2 5t 3 6 + 5t 4 24 t 5 24 375 t 7 56 + . . . . size 12{ {}= { {5t rSup { size 8{2} } } over {2} } - { {5t rSup { size 8{3} } } over {6} } + { {5t rSup { size 8{4} } } over {"24"} } - { {t rSup { size 8{5} } } over {"24"} } - { {"375"t rSup { size 8{7} } } over {"56"} } + "." "." "." "." } {}

Giới hạn chuổi sau số hạn bậc bốn là:

i = 5t 2 2 5t 3 6 + 5t 4 24 i = 5t 2 2 5t 3 6 + 5t 4 24 size 12{i= { {5t rSup { size 8{2} } } over {2} } - { {5t rSup { size 8{3} } } over {6} } + { {5t rSup { size 8{4} } } over {"24"} } } {}

Nếu hàm dùng xấp xỉ i chính xác bốn số thập phân với số hạn xấp xỉ đầu tiên không chú ý đến sai số lớn thì .

5log t  log0,00120

log t  9,415836 - 10

t  0,2605

Giá trị giới hạn là hàm xấp xỉ hợp lý. Vì vậy, trong ví dụ này hàm có thể dùng chỉ để thu được y cho trong khoảng 0  t  0,2; Bởi vì cho t > 0,2 thì e(t) = 1. Cho nên, hàm xấp xỉ khác phải chính xác cho trong khoảng 0,2  t 0,3 như sau:

i = 0, 09367 + 0,2 t 1 i 3i 3 dt i = 0, 09367 + 0,2 t 1 i 3i 3 dt size 12{i=0,"09367"+ Int rSub { size 8{`0,2} } rSup { size 8{`t} } { left (`1 - i - 3i rSup { size 8{3} } right )`} ital "dt"} {}

i ( 1 ) = 0, 09367 + 0,2 t 1 0, 09367 3 0, 09367 3 dt = 0,09367 + 0,90386 ( t - 0,2 ) i ( 1 ) = 0, 09367 + 0,2 t 1 0, 09367 3 0, 09367 3 dt = 0,09367 + 0,90386 ( t - 0,2 ) size 12{i rSup { size 8{ \( 1 \) } } =0,"09367"+ Int rSub { size 8{`0,2} } rSup { size 8{`t} } { left lbrace `1 - 0,"09367" - 3 left (0,"09367" right ) rSup { size 8{3} } right rbrace } ` ital "dt"`=" 0,09367"+"0,90386" \( "t - 0,2" \) } {}

i ( 2 ) = 0, 09367 + 0,2 t 1 0, 09367 0, 90386 t 0,2 3 0, 09367 + 0, 90386 ( t 0,2 ) 3 dt i ( 2 ) = 0, 09367 + 0,2 t 1 0, 09367 0, 90386 t 0,2 3 0, 09367 + 0, 90386 ( t 0,2 ) 3 dt size 12{i rSup { size 8{ \( 2 \) } } =0,"09367"+ Int rSub { size 8{`0,2} } rSup { size 8{``t} } { left lbrace `1 - 0,"09367" - 0,"90386" left (t - 0,2 right ) - 3 left [0,"09367"+0,"90386" \( t - 0,2 \) right ] rSup { size 8{`3} } right rbrace } ` ital "dt"} {} = 0, 09367 + 0, 90386 0,2 t 1 1, 07897 ( t 0,2 ) 0, 76189 t 0,2 2 2, 45089 ( t 0,2 ) 3 dt = 0, 09367 + 0, 90386 0,2 t 1 1, 07897 ( t 0,2 ) 0, 76189 t 0,2 2 2, 45089 ( t 0,2 ) 3 dt size 12{ {}=0,"09367"+0,"90386" Int rSub { size 8{`0,2} } rSup { size 8{`t} } { left lbrace `1 - 1,"07897" \( t - 0,2 \) - 0,"76189" left (t - 0,2 right ) rSup { size 8{2} } - 2,"45089" \( t - 0,2 \) rSup { size 8{3} } right rbrace } ` ital "dt"} {} = 0, 09367 + 0, 90386 x x ( t 0,2 ) 1, 07897 ( t 0,2 ) 2 2 0, 76189 ( t 0,2 ) 3 3 2, 45089 ( t 0,2 ) 4 4 dt = 0, 09367 + 0, 90386 x x ( t 0,2 ) 1, 07897 ( t 0,2 ) 2 2 0, 76189 ( t 0,2 ) 3 3 2, 45089 ( t 0,2 ) 4 4 dt alignl { stack { size 12{ {}=0,"09367"+0,"90386"`x} {} # size 12{~x` left lbrace \( `t - 0,2 \) - 1,"07897" { { \( t - 0,2 \) rSup { size 8{2} } } over {2} } - 0,"76189" { { \( t - 0,2 \) rSup { size 8{3} } } over {3} } - 2,"45089" { { \( t - 0,2 \) rSup { size 8{4} } } over {4} } right rbrace ` ital "dt"} {} } } {}

Cuối cùng, ta có:

i(3) = 0,09367 + 0,90386(t - 0,2) - 0,48762(t - 0,2)2 -

- 0,05420(t - 0,2)3 - 0,30611(t - 0,2)4 + 0,86646(t - 0,2)5 ....

Chuỗi giới hạn, hàm xấp xỉ là:

i = 0,09367 + 0,90386(t - 0,2) -

- 0,48762(t - 0,2)2 - 0,05420(t - 0,2)3 - 0,30611(t - 0,2)4

Cho i hiệu chỉnh trong bốn số thập phân, ta có:

0,86646(t - 0,2)5  0,00005

(t - 0,2)  0,14198

Hàm hợp lý cho trong khoảng 0,2  t 0,342

Giá trị thu được bằng phương pháp Picard được đưa vào trong bảng 2.5.

2.5. SO SÁNH CÁC PHƯƠNG PHÁP.

Trong bài giải của phương trình vi phân hàm quan hệ giữa biến phụ thuộc y và biến độc lập x cần tìm để thỏa mãn phương trình vi phân. Bài giải trong giải tích là rất khó và có một số vấn đề không thể tìm được. Phương pháp số dùng để tìm lời giải bằng cách biểu diễn y như một số hàm của biến độc lập x từ mỗi giá trị xấp xỉ của y có thể thu được bằng sự thay thế hoàn toàn hay biểu diễn tương đương quan hệ giữa các giá trị liên tiếp của y xác định cho việc chọn giá trị của x. Phương pháp Picard là phương pháp số kiểu đầu tiên. Phương pháp Euler, Runge-Kutta, và Milne là ví dụ cho kiểu thứ hai.

Khó khăn chủ yếu phát sinh từ phương pháp xấp xỉ y bằng hàm số, như phương pháp Picard, tìm thấy trong lần lặp lại sự tích phân hiện tại phải thực hiện để thu được hàm thỏa mãn. Vì vậy phương pháp này là không thực tế trong hầu hết các trường hợp và ít được dùng.

Bảng 2.5: Giải bằng phương pháp Picard.

Bảng 5
n Thời gian tn Sức điện động en Dòng điện in
0123456789101112 00,0250,0500,0750,1000,1250,1500,1750,2000,2250,2500,2750,300 00,1250,2500,3750,5000,6250,7500,8751,0001,0001,0001,0001,000 00,001550,006150,013720,024190,037490,053540,072290,093670,115960,137640,158680,17910

Các phương pháp theo kiểu thứ hai đòi hỏi phép tính số học đơn giản đo đó thích hợp cho việc giải bằng máy tính số của các phương trình vi phân. Trong trường hợp tổng quát, đơn giản quan hệ đòi hỏi dùng trong một khoảng nhỏ cho các biến độc lập nhưng ngược lại nhiều phương pháp phức tạp có thể dùng trong khoảng tương đối lớn tốn nhiều công sức trong việc chính xác hóa lời giải. Phương pháp Euler là đơn giản nhất, nhưng trừ khi khoảng tính rất nhỏ thì dùng nó cũng không đúng với thực tế. Phương pháp biến đổi Euler cũng sử dụng đơn giản và có thêm thuận lợi kiểm tra hệ thống vốn có trong quá trình thu được để cải thiện sự ước lượng cho y. Phương pháp có sự chính xác giới hạn, vì vậy đòi hỏi dùng khoảng giá trị nhỏ cho biến độc lập. Phương pháp Runge-Kutta đòi hỏi số rất lớn của phép tính số học, nhưng kết quả cũng không chính xác.

Phương pháp dự đoán sửa đổi của Milne là ít khó khăn hơn phương pháp Runge-Kutta và so sánh được độ chính xác của bậc h5. Vì vậy, phương pháp của Milne đòi hỏi có bốn giá trị ban đầu cho biến phụ thuộc phải thu được bằng một số phương pháp khác, hầu như phương pháp biến đổi Euler hay phương pháp Runge-Kutta, là như nhau. Trong sự ứng dụng máy tính cho phương pháp số. Chương trình đòi hỏi bắt đầu lời giải như phương pháp của Milne. Lời giải tiếp tục dùng công thức khác cho dự đoán và sau đó sửa chữa giá trị của y cung cấp quá trình hệ thống cho kiểm tra tốt bằng sửa chữa ước lượng ban đầu. Nếu sự khác nhau giữa dự đoán và giá trị chính xác là đáng kể, khoảng tính có thể được rút gọn lại. Khả năng trong phương pháp của Milne không có hiệu lực trong phương pháp Runge-Kutta.

Bài tập:

2.1. Giải phương trình vi phân.

dy dx = x 2 y dy dx = x 2 y size 12{ { { ital "dy"} over { ital "dx"} } =x rSup { size 8{2} } - y} {}

Cho 0  t  0,3; với khoảng phương trình 0,05 và giá trị ban đầu x0 = 0 và y0 = 1, bằng các phương pháp số sau đây.

  1. Euler
  2. Biến đổi Euler.
  3. Picard
  4. Xấp xỉ bậc bốn Runge-Kutta
  5. Milne dùng giá trị bắt đầu thu được phương pháp Runge-Kutta

2.2. Giải bằng phương pháp biến đổi Euler hệ phương trình vi phân.

dx dt = 2y dx dt = 2y size 12{ { { ital "dx"} over { ital "dt"} } =2y} {}

dy dt = x 2 dy dt = x 2 size 12{ { { ital "dy"} over { ital "dt"} } = - { {x} over {2} } } {}

Cho 0  t  1,0; Với khoảng phương trình 0,2 và giá trị ban đầu i0 = 0,x0 = 0 và

y0 = 1

2.3. Giải bằng xấp xỉ bậc bốn Runge-Kutta phương trình vi phân bậc hai.

y’’ = y + xy’

Cho 0  x  0,4; Với khoảng phương trình 0,1 và giá trị ban đầux0 = 0,y0 = 1, và y’0 = 0

Collection Navigation

Content actions

Download:

Collection as:

EPUB (?)

What is an EPUB file?

EPUB is an electronic book format that can be read on a variety of mobile devices.

Downloading to a reading device

For detailed instructions on how to download this content's EPUB to your specific device, click the "(?)" link.

| More downloads ...

Module as:

PDF | More downloads ...

Add:

Collection to:

My Favorites (?)

'My Favorites' is a special kind of lens which you can use to bookmark modules and collections. 'My Favorites' can only be seen by you, and collections saved in 'My Favorites' can remember the last module you were on. You need an account to use 'My Favorites'.

| A lens I own (?)

Definition of a lens

Lenses

A lens is a custom view of the content in the repository. You can think of it as a fancy kind of list that will let you see content through the eyes of organizations and people you trust.

What is in a lens?

Lens makers point to materials (modules and collections), creating a guide that includes their own comments and descriptive tags about the content.

Who can create a lens?

Any individual member, a community, or a respected organization.

What are tags? tag icon

Tags are descriptors added by lens makers to help label content, attaching a vocabulary that is meaningful in the context of the lens.

| External bookmarks

Module to:

My Favorites (?)

'My Favorites' is a special kind of lens which you can use to bookmark modules and collections. 'My Favorites' can only be seen by you, and collections saved in 'My Favorites' can remember the last module you were on. You need an account to use 'My Favorites'.

| A lens I own (?)

Definition of a lens

Lenses

A lens is a custom view of the content in the repository. You can think of it as a fancy kind of list that will let you see content through the eyes of organizations and people you trust.

What is in a lens?

Lens makers point to materials (modules and collections), creating a guide that includes their own comments and descriptive tags about the content.

Who can create a lens?

Any individual member, a community, or a respected organization.

What are tags? tag icon

Tags are descriptors added by lens makers to help label content, attaching a vocabulary that is meaningful in the context of the lens.

| External bookmarks