非线性方程组计算
- 格式:doc
- 大小:27.00 KB
- 文档页数:3
在科学与工程计算中,经常遇到求解非线性方程组的问题;非线性方程组在收敛速度及收敛性比线性方程组要差,特别对于非凸的非线性方程组,其求解更是困难。
下面简要介绍非线性方程组的三种解法——牛顿法、拟牛顿法、同伦算法,分析三种解法的适用性,并附Matlab 原程序。
(一)、牛顿迭代法
迭代公式为:x k+1=x k-f(x k)/f'(x k);牛顿迭代法是解非线性方程组比较经典的方法,在局部收敛点附近是平方收敛的;但其解依赖于初始解,且迭代每一步都要计算f'(x k),不仅计算量大而且有时会发生计算困难。
(二)、拟牛顿迭代法
拟牛顿法是为了解决求Jacobi矩阵时带来的困难,现已成为解决非线性方程组和最优化问题的最有效方法之一。
其迭代格式为:
x(k+1)=x(k)-A k-1F(x(k))
A k+1=A k+[(y k-A k s k)(y k-A k s k)T]/[(y k-A k s k)T s k]
在一定条件下,计算H的序列是超收敛的,但稳定性较差,有时迭代效果不理想。
(三)、同伦算法
同伦算法基本思想是从容易求解的方程组开始,逐步过渡到原方程组的求解,从而得到问题的解。
非线性方程组为:F(x)=0,其解为X*。
构造泛函 G:[0,1]XR n->R n
G定义为:G(λ,x)=λ F(x)+(1-λ)[F(x)-F(x(0))]=F(x)+(λ-1)F(x(0))
(其中:x(0)为任意给的初值,假定为λ函数(λ=0))
对于λ的方程G(λ,x)=0,
当λ=0时,0=G(0,x)=F(x)-F(x(0));x(0)是方程的解;
当λ=1时,0=G(1,x)=F(x);x*是方程的解,即x(1)=x*
基于这个思想我们最后可以得到如下关系式:
x'(λ)=-[J(x(λ))]-1F(x(0)) ( 0<=λ<=1,对初始值x(0) )
J为雅可比矩阵,由上面的式子,对λ在[0,1]上积分,就可得到x*=x(1)
上面的非线性方程组问题就转化为数值积分问题。
附:同伦算法的Matlab程序
%% N 为积分部数,h 为积分步长,其值视计算精度而定%% x0 为初始值,可任意给定
%% Equfun(x)为非线性方程组表达式(列向量),
%% Jacobfun(x)为雅可比矩阵
N=100;
h=1/N;
x0=[0 0 0]';
x=x0;
%format long
f=Equfun(x);
b=-h*f;
for i=1:N
A=Jacobfun(x);
k1=inv(A)*b;
A=Jacobfun(x+0.5*k1);
k2=inv(A)*b;
A=Jacobfun(x+0.5*k2);
k3=inv(A)*b;
A=Jacobfun(x+0.5*k3);
k4=inv(A)*b;
x=x+(k1+2*k2+2*k3+k4)/6;
end
disp('The Solution is:')
disp('x=');disp(x);
%%%%%%%%%%%%%%%%%%%%%%%%%
%%%下面是非线性方程组及其雅可比行列式%%%%%%%%%%%%%%%%%%%%%%%%%
function y=Equfun(x)
%nonliner functions
y=[3*x(1)-cos(x(2)*x(3))-1/2;
x(1)^2-81*(x(2)+0.1)^2+sin(x(3))+1.06;
exp(-x(1)*x(2))+20*x(3)+(10*pi-3)/3]; %%%%%%%%%%%%%%%%%%%%%%%%%
function y=Jacobfun(x)
%Jacobi function
%J=zeros(length(x));
y(1,1)=3;
y(1,2)=x(3)*sin(x(2)*x(3));
y(1,3)=x(2)*sin(x(2)*x(3));
y(2,1)=2*x(1);
y(2,2)=-162*(x(2)+0.1);
y(2,3)=cos(x(3));
y(3,1)=-x(2)*exp(-x(1)*x(2));
y(3,2)=-x(1)*exp(-x(1)*x(2));
y(3,3)=20;。