CHƯƠNG 5: CÁC PHƯƠNG TRÌNH PHI TUYẾN

Chia sẻ: Hikaru Hikaru | Ngày: | Loại File: PDF | Số trang:70

Thêm vào BST

Báo xấu

478
lượt xem 107
download

Download Vui lòng tải xuống để xem tài liệu đầy đủ

Nếu phương trình đại số hay siêu việt khá phức tạp thì ít khi tìm được nghiệm đúng. Bởi vậy việc tìm nghiệm gần đúng và ước lượng sai số là rất cần thiết. Ta xét phương trình : f(x) = 0 (1) với f(x) là hàm cho trước của biến x. Chúng ta cần tìm giá trị gần đúng của nghiệm của phương trình này. Quá trình giải thường chia làm hai bước: bước sơ bộ và bước kiện toàn nghiệm. Bước giải sơ bộ có 3 nhiệm vụ: vây nghiệm, tách nghiệm và thu hẹp khoảng chứa nghiệm. Vây nghiệm là tìm xem các nghiệm của phương trình có thể nằm trên ...

Chủ đề:

Bình luận(0) Đăng nhập để gửi bình luận!

Lưu

Nội dung Text: CHƯƠNG 5: CÁC PHƯƠNG TRÌNH PHI TUYẾN

CHƯƠNG 5: CÁC PHƯƠNG TRÌNH PHI TUYẾN §1. KHÁI NIỆM CHUNG Nếu phương trình đại số hay siêu việt khá phức tạp thì ít khi tìm được nghiệm đúng. Bởi vậy việc tìm nghiệm gần đúng và ước lượng sai số là rất cần thiết. Ta xét phương trình : f(x) = 0 (1) với f(x) là hàm cho trước của biến x. Chúng ta cần tìm giá trị gần đúng của nghiệm của phương trình này. Quá trình giải thường chia làm hai bước: bước sơ bộ và bước kiện toàn nghiệm. Bước giải sơ bộ có 3 nhiệm vụ: vây nghiệm, tách nghiệm và thu hẹp khoảng chứa nghiệm. Vây nghiệm là tìm xem các nghiệm của phương trình có thể nằm trên những đoạn nào của trục x. Tách nghiệm là tìm các khoảng chứa nghiệm sao cho trong mỗi khoảng chỉ có đúng một nghiệm. Thu hẹp khoảng chứa nghiệm là làm cho khoảng chứa nghiệm càng nhỏ càng tốt. Sau bước sơ bộ ta có khoảng chứa nghiệm đủ nhỏ. Để xác định khoảng chứa nghiệm ta có thể dùng phương pháp đồ thị. Ngoài ra ta cũng có thể tìm nghiệm bằng phương pháp tìm tăng dần. Ý tưởng của phương pháp này là nếy f1(x).f2(x)
x1 = a; f1 = feval(func,x1); x2 = a + dx; f2 = feval(func,x2); while f1*f2 > 0.0 if x1 >= b x1 = NaN; x2 = NaN; return end x1 = x2; f1 = f2; x2 = x1 + dx; f2 = feval(func,x2); end Khi phát hiện thấy khoảng chứa nghiệm, hàm trả về giá trị biên của đoạn. Nếu không có nghiệm, x1 = x2 = NaN. Ta gọi rootsearch() nhiều lần để phát hiện hết các đoạn chứa nghiệm. Với ví dụ tìm khoảng chứa nghiệm của hàm f(x) = x3 ‐ 10x2 + 5 ta dùng chương trình ctrootsearch.m clear all, clc f = inline(ʹx^3 ‐ 10*x^2 + 5ʹ); [x1, x2] = rootsearch(f,2,10,.2) Bước kiện toàn nghiệm tìm các nghiệm gần đúng theo yêu cầu đặt ra. Có rất nhiều phương pháp xác định nghiệm của (1). Sau đây chúng ta xét từng phương pháp. §2. PHƯƠNG PHÁP LẶP ĐƠN Giả sử phương trình (1) được đưa về dạng tương đương: x = g(x) (2) từ giá trị xo nào đó gọi là giá trị lặp đầu tiên ta lập dãy xấp xỉ bằng công thức: xn = g(xn‐1) (3) với n = 1,2,.... Hàm g(x) được gọi là hàm lặp. Nếu dãy xn → α khi n →∝ thì ta nói phép lặp (3) hội tụ. Ta có định lí: Xét phương pháp lặp (3), giả sử: ‐ [a, b] là khoảng chứa nghiệm α của phương trình (1) tức là của (2) ‐ mọi xn tính theo (3) đều thuộc [a, b] ‐ g(x) có đạo hàm thoả mãn : 242
g′(x) ≤ q < 1 a < x < b (4) trong đó q là một hằng số thì phương pháp lặp (3) hội tụ Ta có thể minh hoạ phép lặp trên bằng hình vẽ sau. x1 xo xo x1 Ta xây dựng hàm simpiter() để lặp function [x, err, xx] = simpiter(g, x0, tolx, maxiter) % giai pt x = g(x) tu x0 bang cah lap %vao : g, x0 = ham va gia tri dau % tolx = sai so mong muon % maxiter = so lan lap max %ra: x = nghiem % err = sai so |x(k) ‐ x(k ‐ 1)| % xx = cac gia tri trung gian if nargin
end x = xx(k); if k == maxiter fprintf(ʹKhong hoi tu sau %d lan lap\nʹ, maxiter) else fprintf(ʹHoi tu sau %d lan lap\nʹ,k) end Để tính lại ví dụ trên ta dùng chương trình ctsimpiter4_2.m clear all, clc f = inline(ʹ‐0.5*((x ‐ 1).^2 ‐ 3)ʹ); [x, ss, xx] = simpiter(f, 0.5,.00001,200) §3. PHƯƠNG PHÁP CHIA ĐÔI CUNG Giả sử cho phương trình f(x) = 0 với f(x) liên tục trên đoạn [a, b] và f(a).f(b)
tol = eps; fa = feval(f, a); fb = feval(f, b); if fa*fb > 0 error(ʹNghiem khong o trong doan nayʹ); end for k = 1: maxiter xx(k) = (a + b)/2; fx = feval(f, xx(k)); err = (b ‐ a)/2; if abs(fx)
Trong đó y = −f(a) 1 h −f(a)+ f(b) (b−a) Tiếp theo dùng cách đó với đoạn [ a, x1] hay [x1, b] mà hai đầu hàm nhận giá trị trái a x1 ξ x dấu ta được nghiệm gần đúng x2 v.v. b Về mặt hình học, phương pháp này có nghĩa là kẻ dây cung của đường cong f(x) qua hai điểm A[a, f(a)] và B[b, f(b)] hay nói cách khác là tuyến tính hoá hàm f(x) trong đoạn [a, b]. Thật vậy phương trình dây cung AB có dạng: f(a) − f(b) af(b) − bf(a) y = x+ a−b a−b Cho x = x1, y = 0 ta có af(b) − bf(a) x1 = (1) f(b) − f(a) Ta xây dựng hàm chord() để thực hiện thuật toán trên function [x, err, xx] = chord(f, a, b, tolx, maxiter) %giai pt f(x) = 0 bg phuong phap day cung. %vao : f ‐ ham can tim nghiem % a/b ‐ khoang tim nghiem % tolx ‐ sai so mong muon cua nghiem % maxiter lan lap max %ra: x ‐ nghiem % err ‐ sai so % xx ‐ cac gia tri trung gian tolfun = eps; fa = feval(f, a); fb = feval(f, b); if fa*fb > 0 error(ʹNghiem khong o trong doan nay !ʹ); end for k = 1: maxiter xx(k) = (a*fb ‐ b*fa)/(fb ‐ fa); %pt.(1) fx = feval(f, xx(k)); err = min(abs(xx(k) ‐ a), abs(b ‐ xx(k))); 246
if abs(fx)
f′′(x i ) 2 e i +1 = − ei y 2f′(xi ) Minh hoạ hình học của thuật toán Newton ‐ Raphson như hình bên. Thuật toán được tóm lược như sau: a x ‐ cho xo x1 b = xo f(x) ‐ tính ∆x = − f ʹ(x) ‐ cho x = x + ∆x ‐ lặp lại bước 2 và 3 cho đến khi |∆x| ≤ ε Ta xây dựng hàm newtonraphson() để thực hiện thuật toán trên. function [x, fx, xx] = newtonraphson(f, df, x0, tolx, maxiter) %giai pt f(x) = 0 bang pp Newton‐Raphson. %vao: f = ftn to be given as a string ’f’ if defined in an M‐file % df = df(x)/dx (neu khong cho se dung dao ham so.) % x0 = gia tri ban dau % tolx = sai so mong muon % maxiter = so lan lap max %ra: x = nghiem % fx = f(x(last)), xx = cac gia tri trung gian h = 1e‐4; h2 = 2*h; tolf = eps; if nargin == 4 & isnumeric(df) maxiter = tolx; tolx = x0; x0 = df; end xx(1) = x0; fx = feval(f,x0); for k = 1: maxiter if ~isnumeric(df) dfdx = feval(df, xx(k)); %dao ham cua ham else dfdx = (feval(f, xx(k) + h)‐feval(f, xx(k) ‐ h))/h2; %dao ham so 248
end dx = ‐fx/dfdx; xx(k+1) = xx(k) + dx; %pt.(3) fx = feval(f, xx(k + 1)); if abs(fx)
Phương trình (2) chính là biểu thức tổng quát của phép lặp. Hai giá trị đầu tiên x1 và x2 cần để khởi động phép lặp. Quá trình lặp được minh hoạ bằng hình a y y f(x) x2 x1 x0 x1 x0 f(x) x x a b Phép lặp có thể không hội tụ (hình b). Tuy nhiên khi hội tụ, nó hội rất nhanh. Ta xây dựng hàm secant() để thực hiện thuật toán trên. function [x, fx, xx] = secant(f, x0, x1, tolx, maxiter) % giai pt f(x) = 0 bg phuong phap day cung %vao : f ‐ ham can tim nghiem % x0, x1 ‐ gia tri khoi dong phep lap % tolx ‐ sai so mong muon % maxiter ‐ so lan lap max %ra: x = nghiem % fx = f(x(last)), xx = cac gia tri trung gian h = (x1 ‐ x0)/100; h2 = 2*h; tolfun = eps; xx(1) = x0; fx = feval(f, x0); for k = 1: maxiter if k
end dx = ‐fx/dfdx; xx(k + 1) = xx(k) + dx; %pt.(2) fx0 = fx; fx = feval(f, xx(k+1)); if abs(fx)
Trong quá trình này ta tính được f(x1), f(x2) và f(x3). Qua 3 điểm này ta có một đường cong bậc 2 và tìm được nghiệm x của đường cong bậc 2 này. Nếu x nằm trong đoạn [x1, x2] như hình trên thì giá trị này được chấp nhận. Tiếp theo ta tìm nghiệm trong đoạn [x1, x3] hay [x3, x2] tuỳ theo vị trí của x. Đa thức nội suy Lagrange qua 3 điểm là: (f − f2 )(f − f3 ) (f − f1 )(f − f3 ) (f − f1 )(f − f2 ) x(y) = x1+ x2+ x3 (f1 − f2 )(f1 − f3 ) (f2 − f1 )(f2 − f3 ) (f3 − f1 )(f3 − f2 ) Cho y = 0 ta có: f f x (f − f ) + f1f3 x 2 (f3 − f1 ) + f1f2 x 3 (f1 − f2 ) x=− 2 3 1 2 3 (1) (f1 − f2 )(f2 − f3 )(f3 − f1 ) Độ thay đổi của nghiệm là: x (f − f )(f − f + f ) + f2 x1 (f2 − f3 ) + f1x 2 (f3 − f1 ) ∆x = x − x 3 = f3 3 1 2 2 3 1 (2) (f1 − f2 )(f2 − f3 )(f3 − f1 ) Ta xây dựng hàm brent() để thực hiện thuật toán trên function root = brent(f, a, b, tol) % giai pt f(x) = 0 bang thuat toan Brent. % Cu phap: root = brent(f, a, b, tol) % vao: f = ham can tim nghiem % a, b = doan chua nghiem % tol = sai so cho truoc if nargin
end if f1*f2 > 0.0 error(ʹNghiem khong nam trong doan nayʹ) end x3 = 0.5*(a + b); % bat dau lap. for i = 1:30 f3 = feval(f, x3); if abs(f3)
end % cho x = x3 & chon x1, x2 moi sao cho x1
. . . . . . . . . . . . . . . . . . . | xn+1 ‐ xn | ≤ q | xn ‐ xn‐1 | Điều này có nghĩa là dãy xi+1 ‐ xi , một cách gần đúng, là một cấp số nhân. Ta cũng coi rằng dãy xn ‐ y với y là nghiệm đúng của (1), gần đúng như một cấp số nhân có công sai q . Như vậy : x n +1 − y = q < 1 (7) xn − y hay : x n +1 − y = q(x n − y) (8) Tương tự ta có : x n + 2 − y = q(x n +1 − y) (9) Từ (8) và (9) ta có : x − x n +1 q = n+2 (10) x n +1 − x n Thay giá trị của q vừa tính ở (10) vào biểu thức của q ở trên ta có : ( x n − x n +1 ) 2 y = xn − (11) x n − 2x n +1 + x n + 2 Công thức (11) được gọi là công thức ngoại suy Adam. Như vậy theo (11) trước hết ta dùng phương pháp lặp để tính giá trị gần đúng xn+2, xn+1, xn của nghiệm và sau đó theo (11) tìm được nghiệm với sai số nhỏ hơn. Khi phương pháp lặp được kết hợp với phương pháp Aitken ta có phương pháp Steffensen. Bắt đầu lặp từ x0, hai bước lặp được dùng để tính: x1 = f(x0) x2 = f(x1) và sau đó dùng thuật toán Aitken để tinh y0 theo (11). Để tiếp tục lặp ta cho x0=y0 và lặp lại bước tính trước. Ta xây dựng hàm aitstef() để thực hiện hai thuật toán trên function [x, y] = aitstef(g, x0, tolx, maxiter) % phuong phap Aitken ‐ Steffenson % giai pt x = g(x) xstart = x0; f0 = 0; f0old = 1.; count = 0; while ((count .0000001)) count = count + 1; f0old = f0; 255
x1 = feval(g, x0); x2 = feval(g, x1); f0 = x0 ‐ (x1 ‐ x0)*(x1 ‐ x0) / (x2 ‐ 2.*x1 + x0); x0 = x1; end x = f0; fprintf(ʹ\nʹ); fprintf(ʹPhuong phap Aitken‐Steffensonʹ); x0 = xstart; count = 0; f0 = 0; x2 = 1.; while ((count .0000001)) count = count+1; x1 = feval(g, x0); x2 = feval(g, x1); f0 = x0 ‐ (x1 ‐ x0)*(x1 ‐ x0) / (x2 ‐ 2.*x1 + x0); x0 = f0; end y = f0; Để tìm nghiệm của phương trình x = (2 ‐ ex + x2)/3 = g(x) ta dùng chương trình ctaitstef.m clear all, clc f = inline(ʹ(2.‐exp(x)+x.^2)/3ʹ); [x, y] = aitstef(f,0., 1e‐4, 50) §9. PHƯƠNG PHÁP MUELLER Trong phương pháp dây cung khi tìm nghiệm trong đoạn [a, b] ta xấp xỉ hàm bằng một đường thẳng. Tuy nhiên để giảm lượng tính toán và để nghiệm hội tụ nhanh hơn ta có thể dùng phương pháp Muller. Nội dung của phương pháp này là thay hàm trong đoạn [a, b] bằng một đường cong bậc 2 mà ta hoàn toàn có thể tìm nghiệm chính xác của nó. 256
Gọi các điểm đó có hoành độ lần lượt là a = x2, b = x1 và ta chọn thêm một điểm x0 nằm trong đoạn [x2, x1]. Gọi h1 = x1 ‐ x0 x0, f0 h2 = x0 ‐ x2 v = x ‐ x0 x1, f1 f(x0) = f0 f(x1) = f1 f(x2) = f2 f(x) x2, f2 h2 h1 h2 γ= h1 Qua 3 điểm này ta có một đường parabol: y = av2 + bv + c Ta tìm các hệ số a, b, c từ các giá trị đã biết v: v = 0 (x = x 0 ) a(0)2 + b(0) + c = f0 v = h1 (x = x1 ) ah1 + bh1 + c = f1 2 v = − h 2 (x = x 2 ) ah 2 − bh 2 + c = f2 2 Từ đó ta có : γf − f (1 + γ ) + f2 a = 1 02 γh 1 (1 + γ ) f1 − f0 − ah 1 2 b= h1 c = f0 Sau đó ta tìm nghiệm của phương trình av2 + bv + c = 0 và có : 2c n = x0 − b ± b 2 − 4ac Dấu của mẫu số được chọn sao cho nó có giá trị tuyệt đối lớn nhất, nghĩa là khi b > 0, ta lấy dấu +, khi b
%vao ‐ f la ham can tim nghiem % ‐ a, b la doan chua nghiem % ‐ maxiter so lan lap max %ra ‐ p xap xi Muller cua f % ‐ y la gia tri y = f(p) % ‐ err sai so thuc cua nghiem. %khoi gan a,b,x0 va cac gia tri tuong ung cua ham x0 = (a + b)/2.; P = [x0 a b]; Y = feval(f, P); delta = 1e‐4; %tinh cac he so cua pt bac hai for k = 1:maxiter h0 = P(1) ‐ P(3); h1 = P(2) ‐ P(3); e0 = Y(1) ‐ Y(3); e1 = Y(2) ‐ Y(3); c = Y(3); denom = h1*h0^2 ‐ h0*h1^2; a = (e0*h1 ‐ e1*h0)/denom; b = (e1*h0^2 ‐ e0*h1^2)/denom; %neu nghiem phuc if b^2‐4*a*c > 0 disc = sqrt(b^2 ‐ 4*a*c); else disc = 0; end %tim nghiem nho nhat if b
P = Q; Y = feval(f, P); end if abs(p‐P(3))
function x = halley1(f, x, maxiter) %ham dung de tim nghiem cua pt f(x) = 0 %vao: ‐ ham can tim nghiem f % ‐ gia tri dau x0 % ‐ so lan lap cuc dai maxiter %ra ‐ nghiem x % dung dao ham so i = 0; h = 1e‐4; while (i