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


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 ĐỀ (*) mọi hệ phương trình. Trong bài báo này
Việc giải xấp xỉ gần đúng một hệ chúng tôi trình bày phương pháp có tên
phương trình phi tuyến là việc mà nhiều tác “từng bước tốt hơn”: biến việc giải hệ
giả trong và ngoài nước đã từ lâu vẫn quan phương trình phi tuyến thành bài toán tìm
tâm. Nhưng chưa có một phương pháp nào cực tiểu của một hàm có dạng tổng bình
tổng quát để giải mọi hệ phương trình. phương. Phương pháp này dựa vào hướng
Hiện nay phương pháp Newton có thể xem gradient giảm dần của hàm để dần xấp xỉ
là phương pháp tốt nhất, tốc độ hội tụ đến cực tiểu của hàm và đó chính là
nhanh nhất để giải hệ phi tuyến, nhưng vẫn nghiệm địa phương của hệ phương trình
không phải là phương pháp vạn năng cho 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:




(*)
TS, Khoa Cơ bản, Đại học Công nghiệp Tp. Hồ
Chí Minh

91
  f1  x1, x 2 ,..., x n  = 0  f1  x1, x 2 ,..., x n  
  
 f 2  x1, x 2 ,..., x n  = 0  f 2  x1, x 2 ,..., x n  
  F = 0 , trong đó F =
 ....................... 
 - - - - - - - - - - -
 
 f x , x ,..., x = 0
  
 f n  x1, x 2 ,..., x n  
 n 1 2 n
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:
n
g  x1 , x 2 ,.., x n  =  fi2 = f12  x1, x 2 ,.., x n  + f 22  x1, x 2 ,.., x n  +.. + f n2  x1, x 2 ,.., x n 
i=1
vậy bài toán đưa về là bài toán tìm p địa phương:

g  p  = min g  x  =0.
xR n

Để 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 x (0) , sau đó tìm điểm dần tốt hơn x (1) sao cho:

 
x (1) = x (0) - α g x (0) 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 đó để: g x (1)  g x (0) .    
Để xác định  , cần tìm cực tiểu của hàm một biến h    g x (0)   g x (0)    .
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à: x (1)  x (0)   g x (0) . Tổng quát ta có  
công thức truy hồi để tìm nghiệm xấp xỉ:

 
x( i 1)  x( i )   g x( i )
Tính công thức gradient hàm g:




92
 n n
  
g  x  =  2fi fi =  2fi fi,1,fi,2 ,fi,3 ,..,fi,n = 2f1 f1,1, f1,2 ,.., f1,n +
i=1 i=1

    
2f 2 f 2,1, f 2,2 ,.., f 2,n + 2f3 f3,1, f3,2 ,.., f3,n +... + 2f n f n,1, f n,2 ,.., f n,n = 
T
 f1,1 f 2,1 M f n,1   f   f1,1 f1,2 M f1,n   f1 
  1   f 
 f1,2 f 2,2 M f n,2  f 2
  = 2  2,1 2,2
f f M f 2,n 
 2
=2
 M M M M   M   M
 M M M  M
   
f M f n,n   f n  f M f n,n  fn 
 1,n f 2,n  n,1 f n,2
= 2JT F  x1, x 2 ,..., x n 
Trong đó J là Jacobi của hàm F  x1 , x2 , ..., xn  . 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: x   x1 , x2 , ..., xn  , số bước N cần tính, sai số chịu đựng tính
T

toán TOL (tolerant)
while k  N do
1. Set g1  g  x1 , x2 , ..., xn 

z  g  x1 , x2 , ..., xn   2J T F  x1 , x2 , ..., xn 
z0  z 2
2. If z0  0 then OUTPUT(“Zero gradient, may have a minimum”) , STOP
z
3. z  , 1  0 , 3  1 , g3  g  x   3 z 
z0
4. While  g3  g1  do steps 5 and 6
5.  3   3 / 2 , g3  g  x   3 z 
6. If  3  TOL / 2 then OUTPUT(“Improvement, may have a
minimum”) , STOP
7.  2   3 / 2 , g2  g  x   2 z 
g  g1 g  g2 h h
8. h1  2 , h2  3 , h3  2 1
 2  1 3  2  3  1
 h1 
9.  0  0.5   2  ,g0  g  x  0 z 
 h3 
10. Tìm  từ  0 ,  3 để g  g  x   z   min g0 , g3
11. x  x  z
12. If g  g1  TOL then OUTPUT(“successful”) , STOP


93
 13. k=k+1
End

2.3. Các ví dụ và so sánh kết quả
 f  6 x  2cos x x  1  0
 1 1 2 3  

Ví dụ 1: giải hệ: 2
 
 f2  9 x2  x1  sin x3  1.06  0.9  0 ,

  x1 x2
 f3  60 x3  3e  10  3  0
n
  1,1,1 . Với hàm g  
(0) T
lấy nghiệm ban đầu x fi2  f12  f22  f32 ,
i 1
Tính Jacobi:

 6 
2sin x2 x3 x3    
2sin x 2 x 3 x 2
 
 x1 2cos  x3  
J=  9 ,
 1  
 x 2  sin x  1.06
3 x1  sin  x3   1.06 
2

 
 x x  x1 x 2 
  3 x2 e 1 2 3 x1 e 60 
f 
 1
và F   f 2 
 
 f3 
 
 k=1,
 -136.9376876467 
Thế x (0)
để tính g1  g x  
 8163.7523; z= 2J F =  24.4584536238
(0)

T ,

 10759.2205558461


tính chuẩn z0  z 2 =10760.1197537; và z 
g x (0)  
 -0.01272640925754 
  0.00227306518733  .
z0  
 0.99991643234926 
Vậy

   
Với 1  0  g1  g x (0)  1z  g x (0)  8163.7523

 3.07635846249709 
Với  3  1  F x  (0)

  3 z =  11.32373711890548 
 29.51313389084277 

94
 
==> g3  g x (0)   3 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
 3.28250948349129 
(0)
 
F x   2 z =  11.48734095680031  ==>

 59.51632652171600 

 
g2  g x (0)   2 z  3684.9269934065

1  2  3
Bây giờ dựa vào bảng mốc nội suy: , để nội suy ra hàm bậc hai theo  :
g1 g2 g3
P    g1  h1 
 1   h3   1    2   g1  h1  h3    2 
1  0 g1  8163.7523
bằng phương pháp Newton tiến nhờ bảng dữ liệu sau:  2  0.5 g2  3684.9269 ;
 3  1 g3  1008.7160
g  g1 g  g2
ta có: h1  2 =-8957.6507317, h2  3 = -5352.42183523,
 2  1 3  2
h h
h3  2 1 =3605.228896495
 3  1
Thế h1 , h3 vào P   , sau đó ta lấy đạo hàm P   để tìm cực tiểu của nó, tức giải
 h 
P /    0 để tìm được nghiệm:   0  0.5   2  h1  =1.4923137322,
 3
 3.34977423120232 

==> F x   0 z = 
(0)
   
11.14453481099873 ==> g0  g x(0)   0 z =135.4224723

 -0.02878932588947 
ta thấy g0 nhỏ hơn g1 và g3 .
 1.01899179529672 
Vậy ta lấy x (1)
x (0)
  0 z   0.99660787360675 
 -0.49218902305463 

 
và tính g x (1)  135.4224723

 58.072382205939 
k=2, tính lại z= 2J F =  203.772128736778  .............
T
 
 -2.042514338131 
x(2)  x(1)  0 z   0.695645, -0.137993, -0.480816 và g x (2)  10.107836......
T
 
95
  
Và thuật toán cứ thế lặp lại với công thức truy hồi: x( i 1)  x( i )   g x( i ) 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à: x( k )  x( k 1)  105 , khoảng 28

bước (xem bảng sau) mới đạt được điều đó.

k x1( k ) x2( k ) x3( k ) 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.522932 -0.202299 -0.527436 0.033205
0
1
0.522660 -0.202218 -0.529208 0.021829
1
1
0.514463 -0.201270 -0.527910 0.014356
2
1
0.514284 -0.201220 -0.529074 0.009448
3
1
0.508890 -0.200672 -0.528222 0.006219
4
1
0.508772 -0.200640 -0.528988 0.004094
5
1
0.505221 -0.200300 -0.528428 0.002696
6
1
0.505144 -0.200279 -0.528932 0.001775
7
1
0.502806 -0.200060 -0.528564 0.001169
8
1
0.502755 -0.200046 -0.528896 0.000770
9
2
0.501215 -0.199904 -0.528653 0.000507
0


96
 2
0.501181 -0.199895 -0.528872 0.000334
1
2
0.500167 -0.199802 -0.528712 0.000220
2
2
0.500145 -0.199796 -0.528856 0.000145
3
2
0.499477 -0.199735 -0.528751 0.000096
4
2
0.499463 -0.199731 -0.528846 0.000063
5
2
0.499023 -0.199691 -0.528777 0.000041
6
2
0.499013 -0.199688 -0.528839 0.000027
7
2
0.498723 -0.199662 -0.528793 0.000018
8

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:

x1( k ) x2( k ) x3( k )
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:

 1 1 
 f  x  cos x x x  1  0
1 2 3 

Ví dụ 2: Giải hệ  f 2  4 1  x1  x2  0.05 x32  0.15 x3  1  0 , lấy x(0)   0, 0, 0  ;
T

 2 2
 f 3   x1  0.1 x2  0.01 x2  x3  1  0

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 x1( k ) x2( k ) x3( k ) 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

97
 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 cứu. Ở đây chúng tôi dựa vào phương pháp
Phương pháp steepest có tốc độ hội tụ tối ưu để tìm cực trị của hàm phi tuyến đa
chậm hơn nhiều so với phương pháp chính biến cho trường hợp riêng khi hàm phi
thống Newton, nhưng phương pháp tuyến đa biến lại là dạng tổng bình phương
steepest luôn chọn được hướng (gradient của các hàm phi tuyến đa biến. Các kết quả
tốt dần) hội tụ của bài toán, thậm chí có này chỉ là bước đầu của việc tìm nghiệm
những bài toán mà phương pháp Newton địa phương, chưa đề cập đến việc đánh giá
không cho ta kết quả, trong khi đó phương sai số giữa nghiệm xấp xỉ và nghiệm chính
pháp steepest lại cho ta kết quả, tất nhiên xác, cũng như chưa xác định vùng cận để
tốc độ hội tụ của phương pháp này rất chọn giá trị ban đầu x (0) sao cho bài toán
chậm. Trên đây cũng chỉ là một trong hội tụ đến nghiệm địa phương cần tìm,
nhiều phương pháp mà thế giới đã nghiên chúng tôi sẽ hẹn những bài báo sau.




98