BỘ GIÁO DỤC VÀ ĐÀO TẠO
TRƯỜNG ĐẠI HỌC SƯ PHẠM HÀ NỘI 2
ĐẶNG THỊ HƯƠNG
PHƯƠNG PHÁP RUNGE-KUTTA GIẢI PHƯƠNG
TRÌNH VI PHÂN THƯỜNG
VÀ PHƯƠNG TRÌNH VI PHÂN ĐẠI SỐ
LUẬN VĂN THẠC SĨ TOÁN HỌC
Chuyên ngành: Toán Giải Tích
Mã số : 60 46 01 02
Người hướng dẫn khoa học
PGS. TS. Tạ Duy Phượng
Hà Nội, 2014
Lời cảm ơn
Luận văn được hoàn thành dưới sự hướng dẫn nhiệt tình và nghiêm túc
của PGS. TS. Tạ Duy Phượng, tôi xin được bày tỏ lòng biết ơn tới Thầy,
đã không chỉ hướng dẫn khoa học, mà còn động viên và khích lệ tác giả
say mê học tập và nghiên cứu.
Tôi cũng xin bày tỏ lòng biết ơn chân thành tới phòng sau đại học, các
thầy cô giáo dạy cao học chuyên ngành Toán giải tích, trường Đại học Sư
phạm Hà Nội 2 đã giúp đỡ tôi trong suốt quá trình học tập.
Tôi cũng xin bày tỏ lòng biết ơn tới Ban giám hiệu, các cán bộ, giáo viên
Khoa cơ bản-Trường Sĩ quan Tăng thiết giáp đã tạo mọi điều kiện thuận
lợi tốt nhất để tôi có thể hoàn thành tốt khóa học cao học.
Cuối cùng tôi cũng xin được gửi lời cảm ơn chân thành tới gia đình, bạn
bè, các đồng nghiệp đã luôn động viên, cổ vũ, giúp đỡ và tạo mọi điều kiện
thuận lợi cho tôi trong quá trình học tập và hoàn thành luận văn.
Hà Nội, tháng 06 năm 2014
Tác giả
Đặng Thị Hương
Lời cam đoan
Tôi xin cam đoan, dưới sự hướng dẫn của PGS. TS. Tạ Duy Phượng, luận
trình vi phân đại số.
Với mong muốn tìm hiểu sâu hơn về phương pháp Runge-Kutta giải phương
trình vi phân, nhằm bổ sung và nâng cao kiến thức đã học trong đại học
và cao học, đồng thời hy vọng có thể sử dụng các kiến thức về giải số
1
2
phương trình trong thực tiễn giảng dạy, tôi chọn đề tài Phương pháp
Runge-Kutta giải phương trình vi phân thường và phương trình
vi phân đại số làm luận văn cao học.
Luận văn gồm hai chương và một phụ lục.
Chương 1 trình bày tổng quan về phương pháp Runge-Kutta giải phương
trình vi phân thường. Lần lượt các phương pháp Runge-Kutta cấp 2, 3, 4
được trình bày, tuy nhiên các bài toán thực tiễn thường đòi hỏi tính toán
theo phương pháp Runge-Kutta cấp cao hơn, nên phương pháp Runge-
Kutta đến cấp 8 được trình bày ở đây.
Chương 2 trình bày phương pháp Runge-Kutta giải phương trình vi phân
ẩn F (t, x, x
) = 0 hoặc phương trình vi phân đại số dạng
x
= f(t, x, y)
0 = g(t, x, y)
theo chương 9, 10 trang 231-299 của [4], chương 5 trang 345 của [8] và một
số tài liệu khác.
Vì giải phương trình vi phân
x
= f(t, x), x ∈ R, t ∈ R, (1.1.1)
với điều kiện ban đầu
x(t
0
) = x
0
, (1.1.2)
tương đương với việc giải phương trình tích phân
x(t) = x
0
+
t
t
0
f(s, x(s))ds, (1.1.3)
nên ta có thể sử dụng quy tắc cầu phương cơ bản để giải phương trình vi
phân.
Mục này trình bày nhất quán cách tìm ra các công thức giải số phương
trình vi phân nhờ quy tắc cầu phương cơ bản.
Nội dung cơ bản của quy tắc cầu phương là: để tính tích phân
b
a
f(t)dt
ta thay f(t) bởi một đa thức nội suy (interpolating polynomial). Tích
phân của hàm f(t) được xấp xỉ bởi tích phân của hàm đa thức (tính được
s
i=1,i=j
(t − c
i
)
c
j
− c
i
.
Khi đó,
b
a
f(t)dt ≈
s
j=1
w
j
f(c
j
),
các trọng số w
j
được tính theo công thức w
j
=
b
t
f(s, x(s))ds. (1.1.5)
Kết hợp với (1.1.4) ta đi đến công thức:
x(t + h) −x(t) = hf(t, x(t)). (1.1.6)
Chia đoạn [a, b] bởi các điểm cách đều a = t
0
< t
1
< . . . < t
n
= b,
gọi h là độ dài bước (stepsize) của biến độc lập t, h =
b − a
n
. Từ công
thức (1.1.6) với kí hiệu x(t
i
) = x
i
, ta có
x
1
= x
0
+ hf(t
0
, x
0
);
x
Suy ra công thức Euler lùi
x
n+1
= x
n
+ hf(t
n+1
, x
n+1
).
Hình 1.2:
Hai phương pháp Euler tiến và Euler lùi là những phương pháp Runge-
Kutta bậc nhất.
• Nếu ta chọn s = 1 và c
1
=
a + b
2
thì ta có công thức xấp xỉ tích phân
bởi diện tích hình chữ nhật ABMN (Hình 1.3):
t+h
t
f(x(s), s)ds = h.f(t +
h
2
, x(t +
h
2
)).
= b thì
L
1
(t) =
t − b
a − b
và L
2
(t) =
t − a
b − a
.
Suy ra w
1
= w
2
=
b − a
2
. Chứng tỏ
b
a
f(t)dt ≈
b − a
2
[f(a) + f(b)].
Như vậy, nếu xấp xỉ tích phân
t+h
pháp ẩn, tức là để tính được giá trị x
n+1
, ta phải giải một phương
trình phi tuyến.
• Nếu chọn s = 3 và c
1
= a, c
2
=
a + b
2
, c
3
= b thì đặt h = b − a, khi đó
ta có:
L
1
(t) =
(t −
a + b
2
)(t − b)
(a −
a + b
2
)(a − b)
=
2
h
2
)(b − a)
=
2
h
2
(t − a)(t −
a + b
2
).
10
Suy ra
w
1
=
b
a
L
1
(t)dt =
b
a
2
h
2
(t −
a + b
2
)(t − b)dt
6
.
Do tính chất đối xứng (hoặc tính trực tiếp), ta có w
3
= w
1
=
h
6
. Từ
các tính toán trên ta đi tới công thức Simpson:
b
a
f(t)dt ≈
h
6
[f(a) + 4f(
a + b
2
) + f(b)].
Suy ra công thức xấp xỉ nghiệm của phương trình vi phân là
Hình 1.5:
x(t + h) −x(t) =
h
6
[f(t, x(t)) + 4f(t +
h
2
, x(t +
)].
Đây là công thức ẩn của phương pháp Runge-Kutta kinh điển cấp
bốn (classical fourth-order Runge-Kutta method).
1.2. Phương pháp Runge-Kutta
1.2.1. Dẫn tới phương pháp Runge-Kutta
Vì phương pháp ẩn đòi hỏi tại mỗi bước phải giải một phương trình phi
tuyến, điều này là không đơn giản, nên ta cố gắng xây dựng các công thức
Runge-Kutta hiển từ công thức hình thang ẩn, công thức điểm giữa ẩn và
công thức Runge-Kutta kinh điển cấp bốn tương ứng như sau.
• Trong công thức hình thang ẩn
x
n+1
= x
n
+
h
2
[f(t
n
, x
n
) + f(t
n+1
, x
n+1
)]
ta thay giá trị x
n+1
ở vế phải bằng công thức Euler tiến
ˆx
h
2
) theo phương pháp
Euler tiến
ˆx
n+1/2
= x
n
+
h
2
f(t
n
, x
n
)
12
và thay vào công thức của phương pháp điểm giữa ẩn
x
n+1
= x
n
+ h.f(t
n
+
h
2
, x(t
n
+
2
+ 2k
3
+ k
4
), n = 0, 1, 2, . . . ,
trong đó k
1
, k
2
, k
3
, k
4
được tính lần lượt:
k
1
= f(t
n
, x
n
);
k
2
= f(t
n
+
h
2
, x
Xét phương trình vi phân
˙x = f(t, x), x ∈ R
n
, t ∈ [0, 1] (1.2.1)
với điều kiện ban đầu
x(0) = x
0
. (1.2.2)
Nội dung cơ bản của phương pháp Runge-Kutta tổng quát như sau. Chia
đoạn [0, 1] thành một lưới đều
t
i
= i.h, i = 0, 1, 2, . . . , N, k =
1
N
và kí hiệu x
i
là giá trị xấp xỉ x(t
i
).
13
Phương pháp Runge-Kutta cho bài toán (1.2.1)-(1.2.2) có dạng
x
i+1
= x
i
+ hΦ(t
i
, x
i
(s-stage là số các điểm mà tại đó phải tính các giá trị của hàm f). Công
thức trên được gọi là công thức Runge-Kutta ẩn tổng quát s-nấc. Phương
pháp Runge-Kutta tổng quát có thể biểu diễn qua bảng Butcher (Butcher
table) sau:
c
1
a
11
a
12
··· a
1s
c
2
a
21
a
22
··· a
2s
··· ··· ··· ··· ···
c
s
a
s1
a
s2
··· a
ss
b
s1
a
s2
··· 0
b
1
b
2
··· b
s
(1.2.6)
Hai phương pháp Runge-Kutta quan trọng thường hay được sử dụng là
phương pháp Runge-Kutta bậc hai và phương pháp Runge-Kutta bậc bốn.
1.2.3. Phương pháp Runge-Kutta bậc hai tổng quát
Giả thiết rằng ta đã biết giá trị của nghiệm x tại t
n
là x
n
. Phương pháp
Runge-Kutta hiển hai nấc cấp hai tổng quát sử dụng điểm (t
n
, x
n
) để xấp
xỉ giá trị của x tại điểm tiếp theo x
n+1
bằng công thức
x
n+1
= x
).
Khái niệm s−nấc (s-stage) thể hiện số lần tính các giá trị của hàm f (tại
các điểm khác nhau trong công thức Runge-Kutta) là s. Để tìm tất cả các
phương pháp Runge-Kutta cấp hai, ta khai triển Taylor bậc nhất của hàm
f(t, x) theo tại (t
0
, x
0
) :
f(t, x) = f(t
0
, x
0
) +
∂f
∂t
(t
0
, x
0
)(t − t
0
) +
∂f
∂x
(t
0
, x
0
)(x − x
21
f(t
n
, x
n
)) = f(t
n
, x
n
) +
∂f
∂t
(t
n
, x
n
)c
2
h
+
∂f
∂x
(t
n
, x
n
)ha
21
f(t
n
n
) +
∂f
∂t
(t
n
, x
n
)c
2
h +
∂f
∂x
(t
n
, x
n
)ha
21
f(t
n
, x
n
)]
= x
n
+ hb
1
f(t
n
n
, x
n
)f(t
n
, x
n
)
= x
n
+ (b
1
+ b
2
)hf(t
n
, x
n
) + b
2
c
2
h
2
∂f
∂t
(t
n
, x
n
) + hx
(t
n
) + h
2
x
(t
n
)
2!
+ ··· .
Vì
x
(t) = f(t, x)
và
x
(t) = (x
(t))
= (f(t, x(t)))
=
∂f
∂t
(t, x(t)) +
(t
n
, x
n
)+
∂f
∂x
(t
n
, x
n
)f(t
n
, x
n
)]+··· .
(1.2.9)
16
So sánh hai cách khai triển (1.2.8) và (1.2.9) theo h ta thấy các hệ số trong
phương pháp Runge-Kutta cấp hai phải thỏa mãn hệ phương trình:
b
1
+ b
2
= 1, c
2
b
2
=
1
2b
2
.
Nếu chọn b
2
=
1
2
thì b
1
=
1
2
và a
21
= c
2
= 1. Khi ấy ta có một phương
pháp Runge-Kutta cấp hai cho phép tính x
n+1
dựa trên công thức
x
n+1
= x
n
+
1
2
hf(t
n
2
. Khi ấy ta có công thức
x
n+1
= x
n
+ h(b
1
k
1
+ b
2
k
2
) = x
n
+ hf(t
n
+
h
2
, x
n
+
h
2
f(t
n
, x
n
1
k
1
+ hb
2
k
2
+ hb
3
k
3
,
17
trong đó:
k
1
= f(t
n
, x
n
);
k
2
= f(t
n
+ c
2
h, x
n
+ ha
n
+ ha
31
f(t
n
, x
n
) + ha
32
[f(t
n
+ c
2
h, x
n
+ ha
21
f(t
n
, x
n
))]
.
Từ bảng Butcher (1.2.6) biểu diễn hệ số của phương pháp Runge-Kutta
hiển ta có điều kiện
c
2
= a
21
∂f
∂t
(t
n
, x
n
)c
2
h
+
∂f
∂x
(t
n
, x
n
)ha
21
f(t
n
, x
n
). (1.2.10)
Đồng thời áp dụng khai triển Taylor bậc hai cho k
2
k
2
= f(t
n
+ c
, x
n
) +
1
2
c
2
2
h
2
∂
2
f
∂t
2
(t
n
, x
n
)
+ c
2
2
h
2
∂
2
f
∂t∂x
(t
) + O(h
3
)
18
Tiếp tục áp dụng khai triển Taylor bậc hai cho k
3
ta có:
k
3
= f(t
n
+ c
3
h, x
n
+ ha
31
f(t
n
, x
n
) + ha
32
f(t
n
+ c
2
h, x
n
+ ha
2
h, x
n
+ ha
21
f(t
n
, x
n
))
+
1
2
c
2
3
h
2
∂
2
f
∂t
2
(t
n
, x
n
)
+ h
2
.f(t
n
+ c
2
h, x
n
+ ha
21
)f(t
n
, x
n
)
+
1
2
h
2
a
2
31
∂
2
f
∂x
2
(t
n
, x
n
h, x
n
+ ha
21
f(t
n
, x
n
))
+
1
2
h
2
a
2
32
∂
2
f
∂x
2
(t
n
, x
n
)f
2
(t
n
)
+ c
2
h
∂f
∂t
(t
n
, x
n
) + c
2
h
∂f
∂x
(t
n
, x
n
)f(t
n
, x
n
)
+
1
2
h
2
c
1
2
h
2
c
2
c
∂
2
f
∂x
2
(t
n
, x
n
)f
2
(t
n
, x
n
) + O(h
3
)
+ h
2
b
3
[f(t
) +
∂f
∂t
(t
n
, x
n
)c
3
h + ha
32
∂f
∂x
(t
n
, x
n
)f(t
n
+ c
2
h, x
n
+ ha
21
f(t
n
, x
n
))
)f(t
n
, x
n
)
+ h
2
c
3
a
32
∂
2
f
∂t∂x
(t
n
, x
n
)f(t
n
+ c
2
h, x
n
+ ha
21
)f(t
n
, x
32
∂
2
f
∂x
2
(t
n
, x
n
)
.f(t
n
, x
n
)f(t
n
+ c
2
h, x
n
+ ha
21
f(t
n
, x
n
))
+
1
3
)].
19
= x
n
+ hb
1
f(t
n
, x
n
) + hb
2
f(t
n
, x
n
) + h
2
b
2
c
2
∂f
∂t
(t
n
, x
n
)
∂t
2
(t
n
, x
n
) + h
3
b
2
c
2
2
∂
2
f
∂t∂x
(t
n
, x
n
)f(t
n
, x
n
)
+
1
2
h
b
3
c
3
∂f
∂t
(t
n
, x
n
) + h
2
b
3
a
31
∂f
∂x
(t
n
, x
n
)f(t
n
, x
n
)
+ h
2
b
3
∂
2
f
∂t
2
(t
n
, x
n
)
+ h
3
b
3
c
3
a
2
31
∂
2
f
∂t∂x
(t
n
, x
n
)f(t
n
n
))
+
1
2
h
3
b
3
a
2
31
∂
2
f
∂x
2
(t
n
, x
n
)f
2
(t
n
, x
n
)
+ h
3
n
))
+
1
2
h
3
b
3
a
2
32
∂
2
f
∂x
2
(t
n
, x
n
)f
2
(t
n
+ c
2
h, x
n
+ ha
b
2
c
2
∂f
∂t
(t
n
, x
n
)
+ h
2
b
2
c
2
∂f
∂x
(t
n
, x
n
)f(t
n
, x
n
)
+
1
, x
n
)f(t
n
, x
n
)
+
1
2
h
3
b
2
c
2
2
∂
2
f
∂x
2
(t
n
, x
n
)f
2
(t
n
(t
n
, x
n
)f(t
n
, x
n
)
+ h
2
b
3
a
32
[f(t
n
, x
n
) +
∂f
∂t
(t
n
, x
n
)c
2
h +
∂f
2
f
∂t
2
(t
n
, x
n
) + h
3
b
3
c
3
a
31
∂
2
f
∂t∂x
(t
n
, x
n
)f(t
n
, x
n
)
+ h
∂
2
f
∂x
2
(t
n
, x
n
)f
2
(t
n
, x
n
) + h
3
b
3
a
31
a
32
∂
2
f
∂x
2
(t
n
(t
n
, x
n
) + o(h
4
)
=x
n
+ h(b
1
+ b
2
+ b
3
)f(t
n
, x
n
) + h
2
∂f
∂t
(t
n
, x
n
)(b
2
c
2
h
3
(b
2
c
2
2
+ b
3
c
2
3
)
∂
2
f
∂t
2
(t
n
, x
n
)f
2
(t
n
, x
n
)
, x
n
)
+ h
3
(
1
2
b
2
c
2
2
+
1
2
b
3
a
2
31
b
3
a
31
a
32
+
1
2
2
f
∂x
2
(t
n
, x
n
)f(t
n
, x
n
) + h
3
b
3
a
32
c
2
∂
2
f
∂t∂x
(t
n
, x
n
). (1.2.11)
Mặt khác, khai triển Taylor bậc ba của nghiệm x để được biểu diễn khác
với
x
(t
n
) = f(t
n
, x
n
)
x
(t) = (x
(t))
= f(t, x(t))
=
∂f
∂t
(t, x(t)) +
∂f
∂x
(t, x(t))f(t, x(t)),
21
ta có
x
(t) = (x
f
∂t
2
(t, x(t)) + 2
∂
2
f
∂t∂x
(t, x(t))f(t, x(t)) +
∂
2
f
∂x
2
(t, x(t))f
2
(t, x(t))
+
∂
2
f
∂t∂x
(t, x(t)) +
∂
2
f
∂x
2
(t, x(t))f(t, x(t)).
Từ đây ta được
) f (t
n
, x
n
)
+
1
6
h
3
∂
2
f
∂t
2
(t
n
, x
n
) +
1
3
h
3
∂
2
f
∂t∂x
(t
n
h
3
∂
2
f
∂t∂x
(t
n
, x
n
)
+
1
6
h
3
∂
2
f
∂x
2
(t
n
, x
n
) f (t
n
, x
n
) (1.2.12)
a
32
=
1
2
1
2
b
2
c
2
2
+
1
2
b
3
c
2
3
=
1
6
b
2
c
2
2
+ b
3
2
b
3
a
31
a
32
+
1
2
b
3
a
2
32
=
1
6
b
3
a
32
a
21
=
1
6
b
3
a