1
Mục lục
Mục lục 1
0.1 Lý do chọn đề tài . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
0.2 Mục đích yêu cầu . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
0.3 Phương pháp nghiên cứu . . . . . . . . . . . . . . . . . . . . . . . . . 3
0.4 Giới hạn của đề tài . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1 Phương pháp số cho phương trình vi phân cấp 1 4
1.1 Phương pháp Euler để giải phương trình vi phân cấp 1 . . . . . . . . 4
1.2 Phương pháp Euler cải tiến (Phương pháp Heun) . . . . . . . . . . . 7
1.3 Phương pháp Runge - Kutta . . . . . . . . . . . . . . . . . . . . . . . 9
1.4 Sai số và sự điều chỉnh bước nhãy RKF (Rung - Kutta - Fehlberg) . . 11
2 Phương pháp nhiều nút 14
2.1 Phương pháp Adams - Bashforth . . . . . . . . . . . . . . . . . . . . 14
2.2 Các phương pháp Adams - Moulton . . . . . . . . . . . . . . . . . . . 15
3 Phương pháp cho hệ phương trình và phương trình vi phân bậc
cao 17
3.1 Phương pháp Euler cho hệ phương trình . . . . . . . . . . . . . . . . 17
3.2 Các phương pháp Runge - Kutta cho hệ phương trình . . . . . . . . . 19
3.3 Phương pháp Runge - Kutta - Nystr¨om (phương pháp RKN) . . . . . 20
4 Phương pháp cho các phương trình vi phân đạo hàm riêng 22
4.1 Các phương pháp cho các phương trình vi phân kiểu Elip . . . . . . . 23
4.1.1 Các phương trình vi phân cho phương trình Laplace và phương
trình Poisson . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
4.1.2 Bài toán Dirichlet . . . . . . . . . . . . . . . . . . . . . . . . . 25
4.1.3 Phương pháp ADI . . . . . . . . . . . . . . . . . . . . . . . . 28
4.2 Bài toán Neumann và nhiễu loạn. Biên không chính quy . . . . . . . 30
4.3 Phương pháp cho phương trình Parabol . . . . . . . . . . . . . . . . . 34
2
4.4 Phương pháp cho phương trình Hypebol . . . . . . . . . . . . . . . . 40
5 Phần kết luận 42


) = 0, và nó có thể
viết tường minh y

= f(x, y). Thường phương trình vi phân người ta cho thêm một
điều kiện ban đầu về nghiệm thoả mãn. Trong phần này ta sẽ xét phương trình vi
phân cấp 1 và điều kiện có dạng:
y

= f(x, y); y(x
0
) = y
0
(1.1)
giả sữ f là bài toán có duy nhất nghiệm trên đoạn có chứa x
0
.
Ta sẽ bàn về các phương pháp để tính các giá trị số của các nghiệm, ở đây ta chỉ
cần nêu một cách thức cho nghiệm của các phương trình không khả dụng hoặc quá
phức tạp để tính toán.
Các phương pháp mà ta thực hiện để làm điều này là phương pháp gần đúng
theo từng bứơc (step-by-step motheds). Ta bắt đầu từ việc có giá trị y
0
= y(x
0
)
và tiến hành theo các bước, tính các giá trị gần đúng của nghiệm y(x) tại các nút
"mesh points"
x
1


(x) = y(x) + hf(x, y)
5
(vế phải đạt được ở trên là từ phương trình vi phân y

(x) = f(x, y)) và theo các quá
trình lặp. Trong bước đầu tiên ta tính
y
1
= y
0
+ hf(x
0
, y
0
)
với cách lấy gần đúng y(x
1
) = y(x
0
+ h). Trong bước thứ hai ta tính
y
2
= y
1
+ hf(x
1
, y
1
)

trong mỗi bước hay nó làm ngắn sai số bộ phận (local truncation error) tỉ lệ
với h
2
, được viết O(h
2
). Trên một khoảng x cố định mà ta muốn giải một phương
trình số các bước tỉ lệ với 1/h. Sai số bộ phận hay sai số toàn phần ở đây tỉ lệ với
h
2
(1/h) = h
1
. Do nguyên nhân này mà phương pháp Euler được gọi là Phương
pháp cấp 1. Ngoài ra còn có sai số do làm tròn và do các phương pháp khác, nó
gây ảnh hưởng đến độ chính xác của các giá trị y
1
, y
2
, ···
Ví dụ 1.1. Dùng phương pháp Euler để giải bài toán điều kiện đầu, chọn h = 0.2
và tính y
1
, ··· , y
5
y

= x + y, y(0) = 0 (1.4)
6
Giải: f(x, y) = x + y, và ta thấy rằng (1.3) trở thành:
y
n+1

, trong miền từ h đến 2h nó được nhân với 2
2
= 4, nhưng vì
ta chỉ cần một nữa của bước nhãy trước đó, nên nó chỉ nhân với 4/2 = 2. Do đó
hiệu số 2
2
−
2
= 0.040, sai số 
2
của y
2
trong bảng 1.1 (với kết quả thực là 0.052)
và 0.114 của y
4
(chính xác là 0.152).
n x
n
y
n
0.2(x
n
+ y
n
) Giá trị đúng Sai số
0 0.0 0.000 0.000 0.000 0.000
1 0.2 0.000 0.040 0.021 0.021
2 0.4 0.040 0.088 0.092 0.052
3 0.6 0.128 0.146 0.222 0.094
4 0.8 0.274 0.215 0.426 0.152

(x
n
) |= T OL, nên h
n
=

2T OL
| y

(x
n
) |
(1.5)
y

(x) phải khác không trên đoạn J : x
0
≤ x ≤ x
N
. Chọn K là giá trị cực tiểu của
| y

(x
n
) | trên J và giả sử K > 0 . Giá trị cực tiểu | y

(x
n
) | tương ứng với giá trị
cực đại h = H =

cao và chính xác. Nhưng đây là vấn đề thực hành. Nếu ta thế y

= f (x, y(x)) vào
(1.2) ta được:
y(x + h) = y(x) + hf +
1
2
h
2
f

+
1
6
h
3
f

+ ··· (1.7)
ở đây, vì f phụ thuộc vào x nên
f

= f
x
+ f
y
y

= f
x

1
2
h[f(x
n
, y
n
) + f(x
n+1
, y
∗
n+1
)] (1.8b)
Hình 2: Phương pháp Euler cải tiến
8
Phương pháp này có cách giải thích hình học đơn giãn. Thực tế ta có thể trong
khoảng từ x
n
đến x
n
+
1
2
h ta có nghiệm gần đúng y bằng đường thẳng qua điểm
(x
n
, y
n
) với hệ số góc f(x
n
, y

Giải: Với bài toán lúc này.
k
1
= 0.2(x
n
+ y
n
)
k
2
= 0.2(x
n
+ 0.2)y
n
+ 0.2(x
n
+ y
n
)
y
n+1
= y
n
+
0.2
2
(2.2x
n
+ 2.2y
n

0
+ Nh, ở đây, f như
một bài toán có duy nhất nghiệm trên đoạn [x
0
, x
N
].
Nhập vào: giá trị ban đầu x
0
, y
0
, bước đi(bước nhảy) h số bước đi trong
khoảng khảo sát N
Giá trị ra: Giá trị gần đúng y
n+1
của nghiệm y(x
n+1
) tại x
n+1
= x
0
+(n+1)h.
Với n = 0, 1, ··· , N − 1.
For n = 0, 1, ··· , N − 1 do:
x
n+1
= x
n
+ h
k

Tạm dừng
Kết thúc Phương pháp Euler
Bảng 1.2. Phương pháp Euler cải tiến(Phương pháp Heun)
9
n x
n
y
n
0.22(x
n
+ y
n
) + 0.02 Giá trị chính xác Sai số
0 0.0 0.0000 0.0200 0.0000 0.0000
1 0.2 0.0200 0.0684 0.0214 0.0014
2 0.4 0.0884 0.1274 0.0918 0.0034
3 0.6 0.2158 0.1995 0.2221 0.0063
4 0.8 0.4153 0.2874 0.4255 0.0102
5 1.0 0.7027 0.7183 0.0156
Bảng 1.3. Phương pháp Euler cải tiến áp dụng cho (1.4) và sai số
Sai số bộ phận: Sai số bộ phận của phương pháp Euler cải tiến là số hạng h
3
.
Thật vậy, đặt

f
n
= f(x
n
, y(x


f
n
+

f
n+1
và dùng khai triển
Taylor ta được từ (1.8b)
y
n+1
− y
n
≈
1
2
h[

f
n
+

f
n+1
]
=
1
2
h[



f
n

+
1
4
h
3

f
n

+ ···
(1.10)
Trừ (1.9) cho (1.10) ta có sai số
h
3
6

f
n

−
h
3
4

f
n

3
, k
4
và tiếp đó tính giá trị mới y
n+1
. Công
thức có vẻ phức tạp ở dạng đầu tiên, nhưng thực tế nó lại dễ cho một bài toán.
Việc tính tay f(x, y) bằng tay thì rất phức tạp nhưng nó sẽ không thành vấn đề đối
với máy tính. Phương pháp này là hoàn hảo cho máy tính bởi vì không cần đến các
tính toán ban đầu, nó làm sáng tỏ sự tích lũy, và dùng lặp lại các kết quả một cách
trực tiếp. Nó là số ổn định. Chú ý, nếu f chỉ phụ thuộc vào x, ta dùng công thức
tích phân Simpson.
10
Ví dụ 1.3. áp dụng phương pháp Runge - Kutta đối với bài toán (1.4) trong ví dụ
1.1, chọn h = 0.2 như trước và tính năm bước.
Giải: Với bài toán này, ta có f(x, y) = x + y nên ta có
k
1
= 0.2(x
n
+ y
n
), k
2
= 0.2(x
n
+ 0.1 + y
n
+ 0.5k
1

ta tìm
được k
3
= 0.222(x
n
+ y
n
) + 0.022, và cuối cùng thay giá trị này vào k
4
ta được
k
4
= 0.2444(x
n
+ y
n
) + 0 .0444. Nếu ta dùng các số trên, công thức y
n+1
trong bảng
1.4 trở thành
y
n+1
= y
n
+ 0.2214(x
n
+ y
n
) + 0.0214. (1.11)
Tất nhiên, các quá trình thay vào này là không tiêu biểu cho phương pháp

, x
N
].
Nhập vào: giá trị ban đầu x
0
, y
0
, bước đi(bước nhảy) h số bước đi trong
khoảng khảo sát N
Giá trị ra: Giá trị gần đúng y
n+1
của nghiệm y(x
n+1
) tại x
n+1
= x
0
+(n+1)h.
Với n = 0, 1, ··· , N − 1.
For n = 0, 1, ··· , N − 1 do:
k
1
= hf(x
n
, y
n
)
k
2
= hf(x

n
+ h, y
n
+ k
3
)
x
n+1
= x
n
+ h
y
n+1
= y
n
+
1
6
(k
1
+ 2k
2
+ 2k
3
+ k
4
)
Lấy ra: x
n+1
, y

thức (1.11).
x y = e
x
− x −1 Sai số
Phương pháp Euler Phương pháp
Euler cải tiến
Phương pháp
Runge - Kutta
0.2 0.021 403 0.021 0.0014 0.000 003
0.4 0.091 825 0.052 0.0034 0.000 007
0.6 0.222 119 0.094 0.0063 0.000 011
0.8 0.425 541 0.152 0.0102 0.000 020
1.0 0.718 282 0.229 0.0156 0.000 031
Bảng 1.6. So sánh sự chính xác của ba phương pháp xét trong bài toán (1.4) với
h = 0.2.
1.4 Sai số và sự điều chỉnh bước nhãy RKF (Runge
- Kutta - Fehlberg)
ý tưởng này dùng tích phân thích hợp cho Rung-Kutta và các phương pháp
khác. Trong bảng 1.4 cho RK (Runge - Kutta), nếu ta tính mỗi bước với bước nhãy
h và 2h, thì ít lâu sau sai số cho mỗi bước bằng 2
5
= 32 nhân với giá trị ban đầu.
Tuy nhiên, ta chỉ có một nữa bước cho 2h, nên số chính xác là 2
5
/2 = 16. ở đây sai
số  xấp xĩ với bước nhãy h bằng cở 1/15 nhân với độ chênh lệch δ = y −

y của các
giá trị tương ứng với cở bước nhãy lần lượt là h và 2h
 ≈

) đến (x
n+1
, y
n+1
). Sự chênh lệch
của việc tính giá trị y tại x
n+1
cho sai số ước tính được dùng cho điều chỉnh cở bước
nhãy h. E. Fhlberg đưa ra hai công thức RK với công thức này chỉ cần 6 giá trị
hàm cho mỗi bước. Ta trình bày dạng các hàm ở đây vì RKF đã trở nên quá phổ
biến. Chẳng hạn, Maple dùng nó (cho hệ thống các phương trình vi phân).
Phương pháp RK bậc năm của Fhlberg là
y
n+1
= y
n
+ γ
1
k
1
+ ··· + γ
6
k
6
(1.14a)
γ =

16
135
0

∗
5
k
6
(1.15a)
γ
∗
=

25
216
0
1408
2565
2197
4104
−
1
5

(1.15b)
Trong cả hai công thức ta chỉ dùng 6 giá trị hàm khác nhau
k
1
= hf(x
n
, y
n
)
k

k
2
k
4
= hf(x
n
+
12
13
h, y
n
+
1932
2197
k
1
−
7200
2197
k
2
+
7296
2197
k
3
) (1.16)
k
5
= hf(x

−
8
27
k
1
+ 2k
2
−
3544
2565
k
3
+
1859
4104
k
4
−
11
40
k
5
13
Sự sai khác của (1.14a) và (1.15a) cho sai số ước tính

n+1
≈ y
n+1
− y
∗

k
1
= 0.200000000000 k
2
= 0.200062500000 k
3
= 0.200140756867
k
4
= 0.200856926154 k
5
= 0.201006676700 k
6
= 0.200250418651
y
1
= 1.20033467253 y
∗
1
= 1.20033466949 
1
= 0.00000000304
Giá trị chính xác đến 12 số thập phân là y(0.1) = 1.20033467209. Do đó sai số chính
xác của y
1
là −4.4 · 10
−10
, nhỏ hơn trong bảng 1.7 bằng hằng số 200.
Bảng 1.8 tổng kết các kết quả cơ bản của các phương pháp trong chương 1
này. Nó có thể các phương pháp này là các sai số ổn định. Chúng là các phương

y
∗
n
y
n
chính xác của y
n
0 0.0 0.000 000 0.000 000 0
1 0.2 0.021 400 0.021 403 3
2 0.4 0.091 818 0.091 825 7
3 0.6 0.222 107 0.222 119 12
4 0.8 0.425 361 0.425 529 0.425 541 12
5 1.0 0.718 066 0.718 270 0.718 282 12
6 1.2 1.119 855 1.120 106 1.120 117 11
7 1.4 1.654 885 1.655 191 1.655 200 9
8 1.6 2.352 653 2.353 026 2.353 032 6
9 1.8 3.249 190 2.249 646 3.249 647 1
10 2.0 4.388 505 4.389 062 4.389 056 −1
Bảng 1.9. Phương pháp Adams - Moulton áp dụng vào bài toán (2.9); giá trị tiên
đoán được tính bằng (2.7) và các giá trị được hiệu chính bằng (2.8).
14
Chương 2
PHƯƠNG PHÁP NHIỀU NÚT
Phương pháp một bước(một nút) là phương pháp mà trong mỗi bước chúng ta
dùng các kết quả đạt được từ các bước trước một cách có thứ tự. Tất cả các phương
pháp trong chương 1 đều là phương pháp một bước. Trái lại, phương pháp nhiều
nút là phương pháp mà mỗi bước dùng các giá trị từ các bước thứ tự. Lý do mà
chúng ta dùng các bước tương đương là thêm các thông tin có thể làm tăng độ chính
xác của kết quả. Các phương pháp đó sẽ lần được trình bày như ở dưới đây:
2.1 Phương pháp Adams - Bashforth

n
) =

x
n+1
x
n
f(x, y(x))dx.
Ta thay f(x, y(x)) bởi đa thức nội suy (interpolation polynomical) p(x), để ta đi
tính tích phân, điều này cho kết quả xấp xĩ y
n+1
của y(x
n+1
) và y
n
của y(x
n
)
y
n+1
= y
n
+

x
n+1
x
n
p(x)dx (2.2)
Chọn các giá trị khác nhau của p(x) sẽ tạo ra các phương pháp khác nhau. Ta giải

, y
n−2
), f
n−3
= f(x
n−3
, y
n−3
). (2.3)
Điều này sẽ đưa ra các biểu thức tốt. Ta có thể tính được p
3
(x) từ biểu thức Newton.
p
3
(x) = f
n
+ r∇f
n
+
1
2
r(r + 1)∇
2
f
n
+
1
6
r(r + 1)(r + 2)∇
3

x
n+1
x
n
p
3
dx = h

1
0
P
3
dr = h

f
n
+
1
2
∇f
n
+
5
12
∇
2
f
n
+
3

− 3f
n−1
+ 3f
n−2
− f
n−3
.
Thế vào (2.4) và gom các số hạng lại ta thư được biểu thức nhiều bước của Phương
pháp Adams - Bashforth bậc bốn:
y
n+1
= y
n
+
h
24
(55f
n
− 59f
n−1
+ 37f
n−2
− 9f
n−3
) (2.5)
Biểu thức này diễn tả giá trị mới y
n+1
[xấp xĩ nghiệm y của (2.1) tại x
n+1
] trong

n−2
. (Trước đây ta đã có tại các điểm x
n
, x
n−1
, x
n−2
, x
n−3
. Ta
tính đa thức nội suy nhưng lúc này ta thay r = (x − x
n+1
h), ta được
p
3
(x) = f
n+1
+ r∇f
n+1
+
1
2
r(r + 1)∇
2
f
n+1
+
1
6
r(r + 1)(r + 2)∇

2
f
n+1
−
1
24
∇
3
f
n+1

.
16
Ta cũng có biểu thức
y
n+1
= y
n
+

x
n+1
x
n
p
3
(x)dx = y
n
+
h

= y
n
+
h
24
(55f
n
− 59f
n−1
+ 37f
n−2
− 9f
n−3
) (2.7)
Giá trị đúng mới y
n+1
đạt được từ (2.6) với f
n+1
được thay bằng f
∗
n+1
= f(x
n+1
, y
∗
n+1
)
vì vậy đạt được giá trị như sau:
y
n+1

nên ta có thể so sánh các kết quả. Ta tính các giá trị ban đầu y
1
, y
2
, y
3
bằng phương
pháp Runge - Kutta cỗ điển. Rồi trong mỗi bước ta tiên đoán bằng cách dùng (2.7)
và hiệu chỉnh bằng cách dùng (2.8) trước khi tiến hành các bước tiếp theo. Các kết
quả trong bảng 1.9 được biễu diển và so sánh với kết quả chính xác. Ta thấy rằng
tính chính xác được cải tiến, khảo sát chính xác.
Nhận xét về sự so sánh của các phương pháp . Biểu thức Adams - Moulton là chính
xác hơn biểu thức Adams - Bashforth cùng bậc. Phương pháp (2.7); (2.8) là phương
pháp số ổn định, trái lại ngoại trừ sữ dụng (2.7) có thể gây ra sự bất ổn định. Sự
điều chỉnh ở mỗi bước sẽ tương đối đơn giản. Nếu |Corrector −P redictor| > T OL,
dùng phép nội suy để phát sinh kết quả cũ tại một nữa tại bước hiện thời và thử
h/2 như bước mới. Phương pháp (2.7), (2.8) chỉ cần 2 giá trị cho mỗi bước, phương
pháp Runge - Kutta cần 4 giá trị; tuy nhiên phương pháp Runge - Kutta có thể có
số bước lớn hơn hai rất lớn.
17
Chương 3
Phương pháp cho hệ phương trình
và các phương trình bậc cao
Ta xét hệ phương trình vi phân cấp một có điều kiện đầu sau:
y

= f(x, y), y(x
0
) = y
0

, ··· , y
m
)
Giả sữ f là bài toán có duy nhất nghiệm theo x trên đoạn chứa x
0
. Hệ phương trình
(3.1) gồm các phương trình vi phân cấp m
y
(m)
= f(x, y, y

, y

, ··· , y
(m−1)
) (3.2)
và điều kiện đầu tương ứng là y(x
0
) = K
1
, y

(x
0
) = K
2
, ··· , y
(m−1)
= K
m

m−1
= y
m
, y

m
= f(x, y
1
, ··· , y
m
). (3.4)
và điều kiện đầu y
1
(x
0
) = K
1
, y
2
(x
0
) = K
2
, ··· , y
m
(x
0
) = K
m
.

1,n
+ hf
1
(x
n
, y
1,n
, y
2,n
)
y
2,n+1
= y
2,n
+ hf
2
(x
n
, y
1,n
, y
2,n
)
và tương tự cho hệ phương trình nhiều hơn hai phương trình. Dựa vào (3.4) các
phương trình đã cho chuyển đổi thành hệ
y

1
= f
1

+ 0.2(−2y
2,n
− 0.75y
1,n
)
Các điều kiện đầu là y(0) = y
1
(0) = 3, y

(0) = y
2
(0) = −2.5 các tính toán được trình
bày trong bảng 3.1. Như các phương pháp đơn giãn, các kết quả sẽ không đủ chính
xác cho mục đích tính toán. Tất nhiên ví dụ chỉ là để minh hoạ cho phương pháp
còn bài toán này ta có thể giải được nghiệm một cách chính xác
y = y
1
= 2e
−0.5x
+ e
−1.5x
, y

= y
2
= −e
−0.5x
− 1.5e
−1.5x
y

3.2 Các phương pháp Runge - Kutta cho hệ phương
trình
Được dùng cho hệ phương trình (3.1) bằng cách viết dưới dạng vectơ, thành
phần m cho m = 1 ép buộc các biểu thức vô hướng trước. Từ phương pháp Runge
- Kutta cổ điển bậc bốn trong bảng 1.4 ta được
y(x
0
) = y
0
, giá trị đầu (3.6)
và cho n = 0, 1, ··· , N − 1 (N số các bước), ta đạt được các tính chất phụ trợ sau
k
1
= hf(x
n
, y
n
), k
2
= hf(x
n
+
1
2
h, y
n
+
1
2
k

y
n+1
= y
n
+
1
6
(k
1
+ 2k
2
+ 2k
3
+ k
4
) (3.8)
Ví dụ 3.2. Phương pháp Runge - Kutta cho hệ phương trình. Phương
trình Airy. Hàm Airy Ai(x)
Giải bài toán sau
y

= xy, y(0) =
1
3
2/3
Γ(2/3)
= 0.35502805, y

(0) = −
1

ở đây f = [f
1
f
2
]
T
trong (3.1) có các thành phần
f
1
(x, y) = y
2
, f
2
(x, y) = xy
1
Ta viết các thành phần (??). Điều kiện đầu (3.6) là y
1,0
= 0.35508805, y
2,0
=
−0.25881940. Ta có thể viết
k
1
= a, k
2
= b, k
3
= c, k
4
= d

h)(y
1,n
+
1
2
a
1
)

20
c = h

y
2,n
+
1
2
b
2
(x
n
+
1
2
h)(y
1,n
+
1
2
b

h. Số thứ hai
trong b là y = y
n
+
1
2
a và thành phần đầu tiên là y
1
= y
1,n
+
1
2
a
1
. Cùng với,
xy
1
= (x
n
+
1
2
h)(y
1,n
+
1
2
a
1

× sai số của y
1
y
2,n
(x
n
)
0 0.0 0.35502 805 0.35502 805 0 -0.25881 940
1 0.2 0.30370 303 0.30370 315 12 -0.25240 464
2 0.4 0.25474 211 0.25474 235 24 -0.23583 073
3 0.6 0.20979 973 0.20980 006 33 -0.21279 185
4 0.8 0.16984 596 0.16984 632 36 -0.18641 171
5 1.0 0.13529 207 0.13529 242 35 -0.15914 687
Bảng 3.2. Phương pháp Runge - Kutta cho hệ phương trình: Các giá trị y
1,n
(x
n
)
của hàm Airy Ai(x) trong ví dụ 3.2
3.3 Phương pháp Runge - Kutta - Nystr¨om (phương
pháp RKN)
Phương pháp RKN là mở rộng trực tiếp phương pháp RK(phương pháp Runge
- Kutta) cho phương trình vi phân cấp hai y

= f(x, y, y

) như đã cho bằng phép
toán E. J. Nystr¨om. Đa thức tốt nhất được biết, ở đây n = 0, 1, ··· , N − 1 (N là
số bước nhảy). Các biểu thức dưới với K =
1

=
1
2
hf(x
n
+
1
2
h, y
n
+ K, y

n
+ k
1
)
k
3
=
1
2
hf(x
n
+
1
2
h, y
n
+ K, y


n+1
= y
n
+ h[y

n
+
1
3
(k
1
+ k
2
+ k
3
)] (3.12)
và xấp xĩ y

n+1
của đạo hàm y

(x
n+1
) trong bước tiếp theo
y

n+1
= y

n

n
)
k
2
= k
3
=
1
2
hf(x
n
+
1
2
h, y
n
+
1
2
h(y

n
+
1
2
k
1
)) (3.14)
k
4

n+1
=
y

n
+
1
2
(
k
1
+ 4
k
2
+
k
4
)
Ví dụ 3.3. Phương pháp Runge - Kutta - Nystr¨om. Phương trình Airy.
Hàm Airy Ai(x)
Với bài toán trong ví dụ 3.2 và h = 0.2 như trước ta đạt được từ (3.14) đơn giản
như sau k
1
= 0.1x
n
y
n
và k
2
= k


n
y(x) chính xác 10
8
× sai số của y
n
0.0 0.355 028 05 -0.258 819 40 0.355 028 05 0
0.2 0.303 703 04 -0.252 404 64 0.303 703 15 11
0.4 0.254 742 11 -0.235 830 70 0.254 742 35 24
0.6 0.209 799 74 -0.212 791 72 0.209 800 06 32
0.8 0.169 845 99 -0.186 411 34 0.169 846 32 33
1.0 0.135 292 18 -0.159 146 09 0.135 292 42 44
Bảng 3.3. Phương pháp Runge - Kutta - Nystr¨om áp dụng cho phương trình Airy,
tính toán của hàm Airy y = Ai(x)
22
Chương 4
PHƯƠNG PHÁP CHO CÁC
PHƯƠNG TRÌNH VI PHÂN
ĐẠO HÀM RIÊNG
Các phương trình vi phân đạo hàm riêng đặc biệt là phương trình Laplace,
phương trình Poisson, phương trình truyền nhiệt và phương trình sóng có cơ sở là
sự phụ thuộc theo thời gian theo các kiểu khác nhau. Các kiểu đó được phân thành
phương trình kiểu elip, phương trình kiểu hypebol và phương trình kiểu parabol.
Chúng được định nghĩa như dưới đây.
Phương trình vi phân đạo hàm riêng được gọi là chuẩn tuyến tính nếu nó tuyến
tính trong đạo hàm cao nhất. Với phương trình chuẩn tuyến tính bậc hai với hai
biến độc lập x, y có thể được viết
au
xx
+ 2bu

=
∂u
∂n
(đạo hàm riêng của u) được quy định trên
C và bài toán thứ ba hay bài toán hỗn hợp nếu u được đặt ra trên một phần của C
và u
n
trên phần còn lại. C là đường cong đóng (hay chứa hai hay nhiều đường).
4.1 Các phương pháp cho các phương trình vi
phân kiểu Elip
4.1.1 Các phương trình vi phân cho phương trình Laplace
và phương trình Poisson
Phương trình Laplace là
∇
2
u = u
xx
+ u
yy
= 0 (4.2)
Phương trình Poisson
∇
2
u = u
xx
+ u
yy
= f(x, y) (4.3)
Các phương trình này là các phương trình elip quan trọng nhất trong áp dụng.
Để đạt được các phương pháp của lời giải số, ta thay đạo hàm riêng tương ứng bằng

h
3
u
xxx
(x, y) + ··· (4.4b)
Lấy (4.4a) trừ cho (4.4b), bỏ các số hạng chứa h
3
, h
4
, |cdots và tính u
x
thì ta được
u
x
(x, y) ≈
1
2h
[u(x + h, y) −u(x − h, y)]. (4.5)
Tương tự,
u(x, y + k) = u(x, y) + ku
y
(x, y) +
1
2
k
2
u
yy
(x, y) + ···
u(x, y −k) = u(x, y) −ku

2
u
xx
(x, y).
Ta tính được
u
xx
(x, y) ≈
1
h
2
[u(x + h, y) −2u(x, y) + u(x −h, y)]. (4.7)
Tương tự
u
yy
(x, y) =
1
k
2
[u(x, y + k) −u(x, y − k)]. (4.8)
u
xy
(x, y) ≈
1
4hk
[u(x + h, y + k) − u(x − h, y + k)
− u(x + h, y − k) + u(x − h, y −k)]. (4.9)
Hình 3 diễn tả các điểm (x + h, y), (x − h, y), cdots trong (4.5), (4.6), (4.7), (4.8),
(4.9).
Hình 3: Các điểm trong (4.5), (4.6), (4.7), (4.8), (4.9).

(4.13)
Hình 4: Các điểm trong phương trình (4.10) và (4.11) và các điểm đã chuyển đổi
cho tiện
Như vậy phương trình (4.10) có thể viết thu gọn như sau



1
1 −4 1
1



u = h
2
f(x, y)
Chú ý rằng (4.11) làm nổi bật u tại (x, y) bằng trung bình của các giá trị của u tại
bốn điểm lân cận. Đây là mô phỏng tính chất trung bình của hàm điều hoà.
4.1.2 Bài toán Dirichlet
Trong nghiệm số của bài toán Dirichlet trong miền R, trước tiên ta chọn h và
đưa vào trong R một lưới gồm các đường nằm ngang và các đường thẳng đứng cách
đều cách nhau một khoảng h. Các điểm giao nhau của chúng gọi là các mắt lưới
(mesh point). Xem hình 5. Rồi ta dùng phương trình sai phân xấp xỉ phương trình
vi phân đạo hàm riêng - dạng (4.11) trong trường hợp của phương trình Laplace -
với các giá trị chưa biết của u tại các mắt lưới trong R đến các điểm khác và các

Trích đoạn Phương pháp cho phương trình Parabol

---