Archive

Posts Tagged ‘System of equations’

Phương pháp Newton cho hệ phi tuyến

October 7, 2009 2 comments

Lý thuyết xấp xỉ nghiệm của hệ phi tuyến tổng quát bằng phương pháp Newton có thể xem trong bài báo dưới đây 

http://benisrael.net/Newton-MP.pdf

Ở đây giới thiệu maplet cho bài xấp xỉ nghiệm của hệ phi tuyến bẳng phương pháp Newton hiểu theo nghĩa cổ điển.

Thí dụ giải hệ phương trình x^3-y+1=0, 2x+y+3=0

Trước hết lập một vài procedure cần thiết

> with(LinearAlgebra): with(MTM):
> F:=proc(F1,F2,x,y) local f1,f2;      # Lập vector hàm
f1:= proc (x,y) options operator, arrow; F1 end proc;
f2:= proc (x,y) options operator, arrow; F2 end proc;
return (Vector([f1(x,y),f2(x,y)]));
end proc:

> JacobianMatrix := proc (F1,F2,x,y) local f1,f2;    # tính ma trận Jacobian
f1:= proc (x,y) options operator, arrow; F1 end proc;
f2:= proc (x,y) options operator, arrow; F2 end proc;
return (jacobian(Vector([f1(x,y),f2(x,y)]),Vector([x,y])))
end proc:

> V := proc (n::nonnegint,F1,F2,x,y)  #  Tìm nghiệm gần đúng bằng đệ quy
if n = 1 then return Vector([-1, 1]);             
         else return (
V(n-1,F1,F2,x,y)-MatrixVectorMultiply(
MatrixInverse(eval(JacobianMatrix(F1,F2,x,y),[x=V(n-1,F1,F2,x,y)[1],y=V(n-1,F1,F2,x,y)[2]])),
eval(F(F1,F2, x,y),[x=V(n-1,F1,F2,x,y)[1],y=V(n-1,F1,F2,x,y)[2]] ))                     );
end if;
end proc:

> U:=proc(F1,F2,x,y,e) local n;  # Tìm nghiệm gần đúng với sai số là  e
n:=1;
while (VectorNorm(V(n+1,F1,F2,x,y)-V(n,F1,F2,x,y))> e) do
           n:=n+1; end do;  
return(evalf(V(n,F1,F2,x,y)));
end proc:

Có thể test một số hệ phương trình hai biến (dĩ nhiên có thể xảy ra trường hợp phương trình vô nghiệm, hoặc không hội tụ, tuy nhiên trong giải thuật ở trên được xin bỏ qua những điều này, việc bổ sung và tổng quát giải thuật không có nhiều khó khăn), chẳng hạn với ví trên

>  U(x^3-y+1, 2*x+y+3, x, y, 1/1000000);

Kết quả:  x=-1.179509057, y=-0.6409818860

Lập một maplet là bước cuối cùng để thu được sản phẩm phần mềm hữu ích
> with(Maplets[Elements]):

> maplet := Maplet(
Window(‘title’ = “Solving System of Equations”, ‘toolbar’ = ‘TB’,
[[“Enter the 1st function with variables x,y”], [TextField[‘f1’](30)],
[“Enter the 2nd function with variables x,y”], [TextField[‘f2’](30)],
[“Enter tolerance”, TextField[‘e’](12)],
[“Solution”],[MathMLViewer[‘ML’]()]]),

ToolBar[TB](
ToolBarButton(caption = “Do It”, image = Image(“run.png”), onclick=Evaluate(‘ML’=’U(f1,f2,x,y,e)’)),
ToolBarSeparator(),
ToolBarButton(caption = “About”, image = Image(“about.png”), onclick = RunDialog(‘DG’)),
ToolBarSeparator(),
ToolBarButton(caption = “Exit”, image = Image(“exit.png”), onclick = Shutdown())),
 
MessageDialog[DG](“Solving System of Equations. Copyringt (C) 2009 Mathematical Diary”)):

> Maplets[Display](maplet);

systemeq

Advertisements