Thuật toán từng bước (có độ dốc) tốt hơn cho việc giải hệ phương trình phi tuyến

TAÏP CHÍ ÑAÏI HOÏC SAØI GOØN Soá 6 - Thaùng 6/2011 
 91 
THUẬT TOÁN TỪNG BƯỚC (CÓ ĐỘ DỐC) TỐT HƠN 
CHO VIỆC GIẢI HỆ PHƯƠNG TRÌNH PHI TUYẾN 
 NGUYỄN PHÚ VINH 
 (*)
TÓM TẮT 
Trong bài báo này chúng tôi trình bày phương pháp có tên “từng bước tốt hơn”: 
biến việc giải hệ phương trình phi tuyến thành bài toán tìm cực tiểu của một hàm có 
dạng tổng bình phương. Phương pháp này dựa vào hướng gradient giảm dần giá trị của 
hàm để dần xấp xỉ đến giá trị cực tiểu của hàm và giá trị đó cũng chính là nghiệm địa 
phương của hệ phương trình ban đầu. Quá trình này đã được lập trình trên Matlab để 
thử nghiệm, so sánh nghiệm và tốc độ hội tụ của phương pháp này với với phương pháp 
Newton. 
ABSTRACT 
In this article, we study the steepest - descent algorithmwhich is a transformation of 
the solving nonlinear equations system into finding the minimum of a multivariable 
function belonging to the form of sum of squares. This method is on the decent of 
gradient of the function. This function value is approached by the minimum of the 
function which is the local solution of the above nonlinear equations system. This 
process is programmed by Matlab in order to test the solution of this problem. It is used 
to compare the solution and the rate of convergence of this method and Newton method. 
1. ĐẶT VẤN ĐỀ (*) 
Việc giải xấp xỉ gần đúng một hệ 
phương trình phi tuyến là việc mà nhiều tác 
giả trong và ngoài nước đã từ lâu vẫn quan 
tâm. Nhưng chưa có một phương pháp nào 
tổng quát để giải mọi hệ phương trình. 
Hiện nay phương pháp Newton có thể xem 
là phương pháp tốt nhất, tốc độ hội tụ 
nhanh nhất để giải hệ phi tuyến, nhưng vẫn 
không phải là phương pháp vạn năng cho 
(*) TS, Khoa Cơ bản, Đại học Công nghiệp Tp. Hồ 
Chí Minh 
mọi hệ phương trình. Trong bài báo này 
chúng tôi trình bày phương pháp có tên 
“từng bước tốt hơn”: biến việc giải hệ 
phương trình phi tuyến thành bài toán tìm 
cực tiểu của một hàm có dạng tổng bình 
phương. Phương pháp này dựa vào hướng 
gradient giảm dần của hàm để dần xấp xỉ 
đến cực tiểu của hàm và đó chính là 
nghiệm địa phương của hệ phương trình 
ban đầu. 
2. PHƯƠNG PHÁP VÀ KẾT QUẢ NGHIÊN CỨU 
2.1. Phương pháp luận về thuật toán töøng böôùc toát hôn (steepest) 
Đặt bài toán: hãy giải hệ phi tuyến đa biến sau: 
 92 
1 1 2 n
2 1 2 n
n 1 2 n
f x , x ,..., x = 0
f x , x ,..., x = 0
F = 0
- - - - - - - - - - -
f x , x ,..., x = 0
, trong đó 
1 1 2 n
2 1 2 n
n 1 2 n
f x , x ,..., x
f x , x ,..., x
F =
.......................
f x , x ,..., x
ta biến bài toán thành bài toán tìm cực tiểu hàm có dạng tổng bình phương sau: 
 1
2 2 2 2
n i 1 1 2 1 2 n n 1 2 n2
n
g x , x ,.., x = f = f x , x ,.., x + f x , x ,.., x +..+ f x , x , .., xn 2
i=1
 
vậy bài toán đưa về là bài toán tìm p địa phương: 
nx R
g p = min g x
=0. 
Để làm giảm từng bước hàm g(x) và dần đến 0 là cực tiểu của hàm, tiến hành với thuật 
toán như sau: bắt đầu với (0)x , sau đó tìm điểm dần tốt hơn (1)x sao cho: 
 (1) (0) (0)x = x -α g x với >0, (chú ý g là hướng có khả năng làm giảm hàm g(x)), 
ta sẽ điều chỉnh như thế nào đó để: (1) (0)g x g x . 
Để xác định , cần tìm cực tiểu của hàm một biến (0) (0)h g x g x  . 
Nhưng việc khảo sát và lấy đạo hàm trực tiếp của h một cách bài bản thì thật là 
khó, thay vào đó làm theo phương pháp steepest như sau: 
Chọn ba số 1 2 3 và hi vọng chúng sẽ gần vị trí cực tiểu của hàm h , bằng 
cách ta xây dựng một tam thức bậc hai P(x) chính là đa thức nội suy Newton tiến của 
 h tại ba điểm mốc nội suy 1 2, và 3 . Sau đó ta tìm  1 3, để P là 
cực tiểu trong đoạn  1 3, và dùng P để xấp xỉ cực tiểu của h . Rồi dùng 
 làm giá trị xấp xỉ mới cho hàm g, tức là: (1) (0) (0)x x g x  . Tổng quát ta có 
công thức truy hồi để tìm nghiệm xấp xỉ: 
 ( 1) ( ) ( )i i ix x g x  
Tính công thức gradient hàm g: 
 93 
   
n n
i i i i,1 i,2 i,3 i,n
i=1 i=1
g x = 2f f = 2f f ,f ,f ,..,f = 1 1,1 1,2 1,n2f f , f ,.., f +
 2 2,1 2,2 2,n2f f , f ,.., f + 3 3,1 3,2 3,n n n,1 n,2 n,n2f f , f ,.., f +...+ 2f f , f ,.., f = 
=
T
1,1 2,1 n,1 1,1 1,2 1,n1 1
1,2 2,2 n,2 2,1 2,2 2,n2 2
n n1,n 2,n n,n n,1 n,2 n,n
f f M f f f M ff f
f f M f f f M ff f
2 = 2
M MM M M M M M M M
f ff f M f f f M f
= T 1 2 n2J F x , x ,..., x 
Trong đó J là Jacobi của hàm 1 2, , ..., nF x x x . Sau đây là thuật toán steepest. 
2.2. Thuật toán từng bước tốt hơn (steepest) 
Nhập giá trị ban đầu: 1 2, , ...,
T
nx x x x , số bước N cần tính, sai số chịu đựng tính 
toán TOL (tolerant) 
while k N do 
 1. Set 1 1 2, , ..., ng g x x x 
 1 2 1 2, , ..., 2 , , ...,
T
n nz g x x x J F x x x  
 0 2
z z 
 2. If 0z 0 then OUTPUT(“Zero gradient, may have a minimum”) , STOP 
 3. 
0
z
z
z
 , 1 0 , 3 1 , 3 3g g x z 
 4. While 3 1g g do steps 5 and 6 
 5. 3 3 / 2 , 3 3g g x z 
 6. If 3 / 2TOL then OUTPUT(“Improvement, may have a 
 minimum”) , STOP 
 7. 2 3 / 2 , 2 2g g x z 
 8. 2 11
2 1
g g
h
, 3 22
3 2
g g
h
, 2 13
3 1
h h
h
 9. 10 2
3
0.5
h
h
, 0 0g g x z 
 10. Tìm từ 0 3, để 0 3min ,g g x z g g 
 11. x x z 
 12. If 1g g TOL then OUTPUT(“successful”) , STOP 
 94 
 13. k=k+1 
End 
2.3. Các ví dụ và so sánh kết quả 
Ví dụ 1: giải hệ: 
1 2
1 1 2 3
2
2 2 1 3
3 3
6 2cos 1 0
9 sin 1.06 0.9 0
60 3 10 3 0
x x
f x x x
f x x x
f x e 
, 
lấy nghiệm ban đầu (0) 1,1,1 Tx . Với hàm 2 2 2 21 2 3
1
n
i
i
g f f f f
  , 
Tính Jacobi: 
 J=
1 2 1 2
2 3 3 2 3 2
3
1
2 2
1 3 1 3
2 1
6 2sin 2sin
2cos
9
sin 1.06 sin 1.06
3 3 60
x x x x
x x x x x x
xx
x x x x
x e x e
, 
 và 
1
2
3
f
F f
f
 k=1, 
Thế (0)x để tính (0)1g g x 8163.7523; z= 2 TJ F =
 -136.9376876467
 24.4584536238
 10759.2205558461
, 
tính chuẩn 0 2
z z =10760.1197537; và 
 (0)
0
 -0.01272640925754
 0.00227306518733
 0.99991643234926
g x
z
z
 
. 
Vậy 
Với (0) (0)1 1 10 g g x z g x 8163.7523 
Với 3 1 F (0) 3x z =
 3.07635846249709
 11.32373711890548
 29.51313389084277
 95 
==> (0)3 3g g x z 1008.71607578 
Vì g3 < g1 nên ta chấp nhận 3 1 và gán 2 0.5 , rồi tính 
F (0) 2x z =
 3.28250948349129
 11.48734095680031
 59.51632652171600
==> 
 (0)2 2g g x z 3684.9269934065 
Bây giờ dựa vào bảng mốc nội suy: 
1 2 3
1 2 3g g g
, để nội suy ra hàm bậc hai theo : 
 1 1 1 3 1 2 1 1 3 2P g h h g h h 
bằng phương pháp Newton tiến nhờ bảng dữ liệu sau: 
1 1
2 2
3 3
0 8163.7523
0.5 3684.9269
1 1008.7160
g
g
g
; 
ta có: 2 11
2 1
g g
h
=-8957.6507317, 3 22
3 2
g g
h
= -5352.42183523, 
 2 13
3 1
h h
h
=3605.228896495 
Thế 1 3,h h vào P , sau đó ta lấy đạo hàm P để tìm cực tiểu của nó, tức giải 
 / 0P để tìm được nghiệm: 10 2
3
0.5
h
h
=1.4923137322, 
==> F (0) 0x z =
 3.34977423120232
 11.14453481099873
 -0.02878932588947
 ==> (0)0 0g g x z =135.4224723 
ta thấy 0g nhỏ hơn 1g và 3g . 
Vậy ta lấy (1) (0) 0
 1.01899179529672
 0.99660787360675
 -0.49218902305463
x x z 
và tính (1)g x 135.4224723 
 k=2, tính lại z= 2 TJ F =
 58.072382205939
 203.772128736778
 -2.042514338131
............. 
 (2) (1) 0 0.695645, -0.137993, -0.480816
T
x x z và (2)g x 10.107836...... 
 96 
Và thuật toán cứ thế lặp lại với công thức truy hồi: ( 1) ( ) ( )i i ix x g x  cho đến 
khi đạt được nghiệm khá gần với nghiệm chính xác là: 
(0.498145, -0.199606, -0.528826)
T 
ở đây ta lấy tiêu chuẩn dừng của thuật toán là: 
( ) ( 1) 5
10
k k
x x
 , khoảng 28 
bước (xem bảng sau) mới đạt được điều đó. 
k ( )
1
k
x 
( )
2
k
x 
( )
3
k
x g(x) 
1 1.018992 0.996608 -0.492189 135.422472 
2 0.695645 -0.137993 -0.480816 10.107836 
3 0.693185 -0.137886 -0.528789 1.814149 
4 0.583604 -0.226612 -0.523367 0.494123 
5 0.582668 -0.225863 -0.530753 0.295721 
6 0.555173 -0.209607 -0.525621 0.182806 
7 0.554550 -0.209322 -0.529834 0.118411 
8 0.535784 -0.204358 -0.526724 0.077100 
9 0.535372 -0.204218 -0.529429 0.050569 
1
0 
0.522932 -0.202299 -0.527436 0.033205 
1
1 
0.522660 -0.202218 -0.529208 0.021829 
1
2 
0.514463 -0.201270 -0.527910 0.014356 
1
3 
0.514284 -0.201220 -0.529074 0.009448 
1
4 
0.508890 -0.200672 -0.528222 0.006219 
1
5 
0.508772 -0.200640 -0.528988 0.004094 
1
6 
0.505221 -0.200300 -0.528428 0.002696 
1
7 
0.505144 -0.200279 -0.528932 0.001775 
1
8 
0.502806 -0.200060 -0.528564 0.001169 
1
9 
0.502755 -0.200046 -0.528896 0.000770 
2
0 
0.501215 -0.199904 -0.528653 0.000507 
 97 
2
1 
0.501181 -0.199895 -0.528872 0.000334 
2
2 
0.500167 -0.199802 -0.528712 0.000220 
2
3 
0.500145 -0.199796 -0.528856 0.000145 
2
4 
0.499477 -0.199735 -0.528751 0.000096 
2
5 
0.499463 -0.199731 -0.528846 0.000063 
2
6 
0.499023 -0.199691 -0.528777 0.000041 
2
7 
0.499013 -0.199688 -0.528839 0.000027 
2
8 
0.498723 -0.199662 -0.528793 0.000018 
Vậy sau 28 bước ta có nghiệm xấp xỉ là: (498145, -0.199606, -0.528826)T 
So với phương pháp Newton chính thống với ví dụ trên như bảng sau: 
( )
1
k
x 
( )
2
k
x 
( )
3
k
x 
1.127638 -0.270927 -0.513022 
0.498513 -0.192263 -0.523877 
0.498150 -0.199606 -0.528826 
0.498145 -0.199606 -0.528826 
Ta thấy phương pháp Newton chính thống có tốc độ hội tụ nhanh hơn nhiều so với 
phương pháp steepest. Ta xét thêm ví dụ sau: 
Ví dụ 2: Giải hệ 
 1 1 1 2 3
24
2 1 2 3 3
2 2
3 1 2 2 3
cos 1 0
1 0.05 0.15 1 0
0.1 0.01 1 0
f x x x x
f x x x x
f x x x x
, lấy (0) 0,0,0 Tx ; 
Kết quả được chương trình in ra như bảng sau: chỉ cần năm bước là hội tụ 
k ( )
1
k
x 
( )
2
k
x 
( )
3
k
x g(x) 
1 0.000000 0.009944 0.994385 0.008090 
2 -0.020018 0.090056 0.995265 0.000449 
3 -0.000128 0.094959 1.000282 0.000025 
 98 
4 -0.001133 0.099438 0.999764 0.000001 
5 -0.000003 0.099718 0.999995 0.000000 
Nghiệm đúng của hệ này là: x1=0, x2=0.1, x3=1. 
3. KẾT LUẬN 
Phương pháp steepest có tốc độ hội tụ 
chậm hơn nhiều so với phương pháp chính 
thống Newton, nhưng phương pháp 
steepest luôn chọn được hướng (gradient 
tốt dần) hội tụ của bài toán, thậm chí có 
những bài toán mà phương pháp Newton 
không cho ta kết quả, trong khi đó phương 
pháp steepest lại cho ta kết quả, tất nhiên 
tốc độ hội tụ của phương pháp này rất 
chậm. Trên đây cũng chỉ là một trong 
nhiều phương pháp mà thế giới đã nghiên 
cứu. Ở đây chúng tôi dựa vào phương pháp 
tối ưu để tìm cực trị của hàm phi tuyến đa 
biến cho trường hợp riêng khi hàm phi 
tuyến đa biến lại là dạng tổng bình phương 
của các hàm phi tuyến đa biến. Các kết quả 
này chỉ là bước đầu của việc tìm nghiệm 
địa phương, chưa đề cập đến việc đánh giá 
sai số giữa nghiệm xấp xỉ và nghiệm chính 
xác, cũng như chưa xác định vùng cận để 
chọn giá trị ban đầu (0)x sao cho bài toán 
hội tụ đến nghiệm địa phương cần tìm, 
chúng tôi sẽ hẹn những bài báo sau.