不动点迭代法求解非线性方程组
不动点迭代法求解非线性方程组
摘要：一般非线性方程组可以写成()0F x =的形式，其中:n
m
F R R →是定义在区域
n D R ?上的向量函数。把方程组()0F x =改写成与之等价的形式：(x G x =)。因此，求
方程组()0F x =的解就转化为求函数的(G x )的不动点。本文首先介绍了多变量函数()F x 的微积分性质，接着介绍了用不动点迭代法求解非线性方程组。 关键词：多变量函数；微积分；不动点
Fixed Point Iteration Method For Solving Nonlinear Equations
Abstract ：General nonlinear equations can be written in the form of ()F x θ=, where the vector function :n m F R R →is defined on the region n D R ?. Transform the equations ()0F x = into its equivalent form: (x G x =).Therefore, we can get the solution of
()0F x = by finding the fixed point of (G x ).In this paper, we first introduce some knowledge
about multivariable calculus, then introduce the fixed point iteration method for solving nonlinear equations.
Key words: multi-variable function; calculus; fixed point
1 引言
一般非线性方程组及其向量表示法:含有n 个方程的n 元非线性方程组的一般形式为
11221212(,...,)0(,...,)0......(,...,)0
n n
m n f x x x f x x x f x x x =??=??
?
?=? （1） 其中，()1,2,...,i f i n =是定义在n
D R ?上的n 元实值函数，且i f 中至少有一个是非线性
函数。
令12
...n x x x x ??
?
?= ?
???，()()()12()...m f x f x F x f x ??
?
?
= ? ?
??
?
，则方程组可以表示为()F x θ= （2） 其中:n
m
F R R →是定义在区域n D R ?上的向量函数。若存在*
x D ∈，使
*()F x θ=，则称*x 是方程组（1）或（2）的解。把方程组()0F x =改写成与之等价的形
式：(x G x =)
其中n n G D R R ?→：。若*x D ∈满足*
*
(x G x =)，则称*
x 为函数(G x )的不动点。因此
(G x )的不动点就是方程组()0F x =的解，求方程组()0F x =的解就转化为求函数的(G x )
的不动点。适当选取初始向量(0)
x
D ∈，构成迭代公式(+1)()(),
0,1,2,...k k x G x k ==迭代
公式也称为求解方程组()0F x =的简单迭代法，又称不动点迭代法。(G x )称为迭代函数。 由于F 是多变量函数，所以我们先考虑多变量函数的微积分性质。
2多变量函数的微积分性质
在之前我们已经学习过很多关于单变量函数的微积分的性质，由于解非线性方程组经常用到的是多变量函数的相关性质，因此我们考虑多变量函数的微积分性质。相对于单变量函数的微积分的性质，多变量函数的微积分性质一些是类似的，一些是不同的。相对于单变量函数的可微的定义，我们事先给出多变量函数的可微定义。
函数:,1n f R R n ??
→>的微积分性质 设函数多变量函数:,1n
f R R n ??
→>。我们首先考虑当f 是连续的函数的情况，如果f 关于n 个变量的偏导数都存在并且连续，把这n 个偏导数组成一个n 维向量，则我们把这个n 维向量称作多变量函数f 的梯度。
定义1：连续可微函数:n
f R R ??
→，如果()i
f
x x ??，1,...,i n =；存在并且连续，则称函数f 在点n
x R ∈上连续可微，并且称()()()1,...,T
n f f f x x x x x ?????=??????
为函数f 在点
n x R ∈的梯度。
如果函数f 在开区域n D R ?上每一点连续可微，则称函数f 在开区域n D R ?连续可微，记作()1f C D ∈。
下面我们给出关于多变量函数f 的梯度的一些性质：
引理1 设:n
f R R ??
→在开凸集n D R ?连续可微，则对于x D ∈以及任意一个非零扰动n
p R
∈，则函数f 在点x 在方向p 上的方向导数定义为
()()()0lim
f x p f x f x p εεε
→+-?≡?存在并且等于()T
f x p ?。对于,x x p D ?+∈，()()()1
()()x p
T
x
f x p f x f x tp pdt f x f z dz ++=+?+≡+???
，并且存在
(),z x x p ∈+使得，()()()T
f x p f x f z p +=+?。
下面我们给我这个引理的证明过程，主要思想是把多变量函数转化为单变量函数，然后利用我们已知的单变量函数微积分的性质来证明多变量函数微积分的性质。
证明：首先在点x 到点x p +的连线上对函数f 进行参数化，转变成单变量函数g 。定义
()x t x tp =+，:,()()g R R g t f x tp →=+。由链式法则，对于01α?≤≤，
()()()()()()()()1
n i i i f x t dx t dg
x dt x t dt ααα=?=?∑ ()()1n i i i f x p x α=?=??∑ ()f x p p α=?+。 因为
()()00()
()lim (0)t g t g f x x g p t
→-?'==?，所以令0α=，我们就可以得到
()()T
f x f x p p
?=??。 由单变量函数的牛顿定理我们可知，1
(1)(0)()g g g t dt '=+?
。根据前面对函数g 的定
义，上式也可以写成()()
1
()T
f x p f x f x tp pdt +=+?+?。这就得到我们所要的证明。最
后，由单变量函数的积分中值定理
()()1
(),0,1g t dt g ξξ''=∈?
，根据函数g 的定义，我们
可以写成()()()T
f x p f x f x p p ξ+=+?+，()0,1ξ∈。对()
1
T
f x tp pdt ?+?进行变量
替换z x tp =+，可得
()()1
T
x p
x
f x tp pdt f z dz +?+=??
?
，从而得证。
函数:n f R R ??
→的微积分性质 下面给出多变量函数二次可微的定义，并进一步给出函数f 的Hessian 矩阵的定义。
定义2：连续可微函数:n
f R R ??→，如果()2i j
f
x x x ???，1,i j n ≤≤存在并且连续，
则称函数:n
f R R ??
→在点x 上二次连续可微；定义一个n n ?矩阵，其中第,i j 元素为()22
()ij i j
f
f x x x x ??=??，1,i j n ≤≤，则称这个矩阵为函数f 的Hessian 矩阵。
如果函数f 在开区域n D R ?上每一点连续可微，则称函数f 在开区域n
D R ?连续可微，记作()2
f C
D ∈。
类似的我们给出关于多变量函数f 的二阶连续可微的一个引理。
引理2：设函数:n
f R R ??
→在开凸集n D R ?二次连续可微，则对于x D ∈以及任意一个非零扰动n
p R ∈，则函数f 在点x 在方向p 上的二阶方向导数
()()220()lim f f
x p x f
p p
x p
εεε
→??+-???=?存在，并且等于2
()T p f x p ?。对于对于
,x x p D ?+∈，存在(),z x x p ∈+使得()2
1()()()2
T T f x p f x f x p p f z p +=+?+
?。 定理的证明过程与一阶连续可微情况的证明过程类似。从Hessian 矩阵的定义可知，只要函数f 是二次连续可微的，那么Hessian 矩阵是对称的。
函数:n m F R R →的微积分性质
我们进一步考虑更复杂的情况，也就是从n
R 空间到m
R 空间的函数，设函数
:n m
F R R →，具体可以写成()112112
(,...,)..().
...(,...,)()n m n m f x x x f x F x f x x x f x ???? ? ?
? ?
? ?== ? ?
? ?
?
?????。其中，非线性联立方程问题是m n =的情况；非线性最小二乘问题是m n ≠的情况。下面我们给出函数F 的相关可微性质：
定义3 连续函数:n m F R R →，如果每一个部分函数,1,...,i f i m =在点m
x R ∈连续可微，
则称函数F 在点m
x R ∈连续可微。函数F 在点x 的导数叫作F 在点x 的Jacobian 矩阵，它
的转置叫作F 在点x 的梯度。通常的表示为()m n F x R ?'∈，()()i
ij j
f F x x x ?'=
?，()()()T
F x J x F x '==?。
如果F 在开区域n D R ?上每一点连续可微，则称函数F 在开区域n D R ?连续可微，记作()1F C D ∈。
用下面一个例子来具体说明这个定义。
例1 设2
2
:F R R →，()112x
f x e x =-，()2
2122f x x x =-。
解： ()()()()1
1
1
1
2
2
2
1121()22x f f x x x x e J x f f x x x x x ???? ?????
-
?== ? ???-??
?????
下面我们研究单实值函数和向量值函数不同方面，对于实值函数存在中值定理，而对于向量值函数，中值定理不一定成立。也就是说，不一定存在n
z R ∈，使得
()()()F x p F x J z p +=+。直观上来看，尽管每个函数i f 满足
()()()T
i i i i f x p f x f z p +=+?，但是点i z 是不同的。以上面例子中的函数来考虑， 11111110122z e e z -??-??????=+
? ? ? ?--??????
??，也就是，11z
e e =-并且121z =，这是不可能的，所以 不存在n
z R ∈，使得()1,1(0,0)()(1,1)T
F F J z =+。
尽管标准的中值定理是不可能的，我们给出一个近似的中值定理，主要是利用牛顿定理和线性积分的三角不等式。其中，单变量向量值函数的积分可以理解为对每一个部分函数进行黎曼积分。
引理3：设:n
m
F R R →在开凸集n D R ?上连续可微，对于,x x p D ?+∈，有
10
()()()()x p
x
F x p F x J x tp dt F z dz +'+-=+≡??
。
上式可以写成如下中值定理的形式：
()()11
00()()F x p F x J x tp pdt J x tp dt p ??+-=+=+????
??。
因此，我们主要介绍了三种函数，从R R →的函数、从n
R R →的函数以及从
n m R R →的函数的可微性质。
3不动点迭代法
把方程组()0F x =改写成与之等价的形式：(x G x =)。其中n
n
G D R R ?→：。若
*x D ∈满足**(x G x =)，则称*x 为函数(G x )的不动点。因此(G x )的不动点就是方程组()0F x =的解，求方程组()0F x =的解就转化为求函数的(G x )的不动点。
适当选取初始向量(0)
x
D ∈，构成迭代公式(+1)()(),
0,1,2,...k k x G x k ==迭代公式也
称为求解方程组()0F x =的简单迭代法，又称不动点迭代法。(G x )称为迭代函数。
定理1（压缩映射原理）设n n
G D R R ?→：在闭域0D D ?上满足条件：
（1）G 把0D 映入它自身，即00()G D D ?；
（2）G 在0D 上是压缩映射，即存在常数(0,1)L ∈，使对任意的0,x y D ∈，
()()G x G y L x y -≤-
则以下结论成立： （1）对任取的(0)
0x
D ∈，由迭代公式(+1)
()()k k x
G x =，0,1,2,...k =产生的序列{}
()0k x D ? 收敛于函数()G x 在区域0D 内存在唯一的不动点*
x ；
（2）成立误差估计式*()
(1)(0) 1k k L x x x x L -≤--，*()()(1) 1k k k L x x x x L
--≤--。 证明：由于(0)
0x
D ∈以及条件（1）可知序列{}
()0k x D ?，又由条件（2）可得
(1)()()()()(1)(1)(0)()()...k k k k k k k x x G x G x L x x L x x +--=-≤-≤≤- 当1m ≥时有
()
()
()
(1)
1(1)(0)
1
1
m
m
k m k k i k i k i i i x
x
x
x
L x x +++-+-==-≤-≤-∑∑(1)(0)
(1)(0)
(1)11k m k L L L x x x x L L
-=-≤--- （1）
因为01L <<，所以当k →∞时，上式的最后一项是无穷小量，由Cauchy 收敛原理，序列
{}()k x 在n R 中收敛，又由0D 是闭区域{}()k x 的极限(*)0x D ∈，由条件（2）知，
()G x 在0D 上连续，因而*(1)()*lim lim ()()k k k k x x G x G x +→∞
→∞
===，即*x 是方程组
()0F x =的解。
设**
0,x y D ∈是()0F x =的两个不同的解，则有
********()()x y G x G y L x y x y -=-<-<-，
这表明*
x 是()0F x =在内的唯一解。 让m →∞，得*()
(1)(0) 1k
k L x x
x x L -≤-- ()
()
()
(1)
()(1)1
1
m
m
k m k k i k i i k k i i x
x
x
x
L x x +++--==-≤-≤-∑∑
()(1)()(1)
(1)11m k k k k L L L x x x x L L
---=-≤--- 说明：（1）简单迭代法的精度控制与终止条件
()(1)
()
k k k x x x
ε--≤；
（2）由
*(+1)*
()
k k x x L x x
-≤-知简单迭代法是线性收敛的；
（3）对线性方程组迭代函数()G x Bx d =+，有=<1L B ；
定理2（局部收敛定理）设*
int()n n G D R R x D ?→∈：,是方程组()0 F x =的解，在*
x
可微。若()G x '谱半径()()1G x ρ'<，则存在开球{
}
*
0,0D x x x
δδ=-<>，对任取的
(0)0x D ∈，由迭代公式(+1)()()k k x G x =，0,1,2...k =产生序列{}
()0k x D ?收敛于*x 。
下面我们给出一个例子，通过来求函数(G x )的不动点来解非线性方程组。 例2 用简单迭代法求解以下方程
112212
3cos sin 0
4sin cos 0x x x x x x --=??
--=? 要求满足精度()(1)
12()
()10k k k x x e k x
--∞
∞
-=
≤
解：设12x x x ??
= ???
，12121(cos sin )3()1(sin cos )4x x G x x x ??
+ ?= ? ?+ ???
则方程组可以改写成()x G x =，并且对于任意的2
2
,x R y R ∈∈，
112212221(cos cos sin sin )3()()1(sin sin cos cos )4x y x y G x G y x y x y ∞
??
-+- ?
-= ? ?-+- ???
()
7
12
A x y A x y
x y ∞∞
∞∞
=-≤-≤
- 其中，1
2121
1sin cos 33
11cos sin 4
4A ξξηη??- ?
=
? ?- ???
，712A ∞∴≤
因此，任取初始向量(0)
x
R ∈，简单迭代法产生序列{}()k x 收敛于原方程组的唯一解。
迭代公式()()(1)()()
112(1)()()212
1=cos sin 3
1sin cos 4
k k k k k k x x x x x x ++?+????=+??
计算结果：
参考文献：
[1] 李庆杨, 关治, 白峰杉. 数值计算原理[M].清华大学出版社.2000;
[2] 李庆杨，王能超，易大义.数值分析[M].武汉：华中理工大学出版社，1986.
[3] 黄象鼎，曾钟钢，马亚男.非线性数值分析的理论与方法[M].武汉：武汉大学出版社，2004.
[4] 曾金平，李湘良.数值计算方法[M].长沙：湖南大学出版社，2004.
[5] 马吕凤，林伟川.现代数值计算方法[M].北京：科学教育出版社，2008.
[6] 王德人.非线性方程组解法与最优化方法。人民教育出版社，1979
[7] 林成森.数值计算方法（第二版）[M].北京：科学出版社，2005.
MATLAB代码 解线性方程组的迭代法
解线性方程组的迭代法 1.rs里查森迭代法求线性方程组Ax=b的解 function[x,n]=rs(A,b,x0,eps,M) if(nargin==3) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值elseif(nargin==4) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1; %迭代过程 while(tol>eps) x=(I-A)*x0+b; n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 2.crs里查森参数迭代法求线性方程组Ax=b的解 function[x,n]=crs(A,b,x0,w,eps,M) if(nargin==4) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值 elseif(nargin==5) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1; %迭代过程 while(tol>eps) x=(I-w*A)*x0+w*b; n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x;
if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 3.grs里查森迭代法求线性方程组Ax=b的解 function[x,n]=grs(A,b,x0,W,eps,M) if(nargin==4) eps=1.0e-6;%eps表示迭代精度 M=10000;%M表示迭代步数的限制值 elseif(nargin==5) M=10000; end I=eye(size(A)); n=0; x=x0; tol=1;%前后两次迭代结果误差 %迭代过程 while(tol>eps) x=(I-W*A)*x0+W*b;%迭代公式 n=n+1;%n为最终求出解时的迭代步数tol=norm(x-x0); x0=x; if(n>=M) disp('Warning:迭代次数太多，可能不收敛！'); return; end end 4.jacobi雅可比迭代法求线性方程组Ax=b的解 function[x,n]=jacobi(A,b,x0,eps,varargin) if nargin==3 eps=1.0e-6; M=200; elseif nargin<3 error return elseif nargin==5 M=varargin{1}; end D=diag(diag(A));%求A的对角矩阵 L=-tril(A,-1);%求A的下三角阵
非线性方程组迭代解法
非线性方程组迭代解法 不动点法（unmovepoints.m） %非线性方程组的不动点法 function [x,n]=unmovepoints(fun,x0,eps) if nargin<3 eps=1e-3; end x1=feval(fun,x0); n=1; while(norm(x1-x0))>=eps x0=x1;x1=feval(fun,x0); n=n+1; if n>100000 disp('无法收敛!'); break end end x=x1; Newton迭代法（newtons.m） %非线性方程组的Newton迭代法 function [x,n]=newtons(fun1,fun2,x0,eps) if nargin<4 eps=1e-3; end x1=x0-feval(fun1,x0)/feval(fun2,x0); n=1; while norm(x1-x0)>=eps x0=x1;x1=x0-feval(fun1,x0)/feval(fun2,x0); n=n+1; if n>100000 disp('无法收敛!'); break end end x=x1;
注：方程组的迭代与方程迭代不同之处在于收敛的判断不能用abs 而应用norm （范数，默认值为向量各元素的平方和的开方；norm(x1-x0)即为向量x1与x0对应元素差的平方和的开方。在对应的函数程序中应注意向量的运算与数量运算的区别。） 用以上方法求解下列非线性方程组： ()0cos 2.0sin 7.02111=--=x x x X f ()0sin 2.0cos 7.02122=+-=x x x X f 函数： %非线性方程组函数(适用于不动点法) function f=nonlinerequs1(x) f(1)=0.7*sin(x(1))+0.2*cos(x(2)); f(2)=0.7*cos(x(1))-0.2*sin(x(2)); %非线性方程组函数(适用于Newton 迭代法) function f=nonlinerequs2(x) f(1)=x(1)-0.7*sin(x(1))-0.2*cos(x(2)); f(2)=x(2)-0.7*cos(x(1))+0.2*sin(x(2)); %非线性方程组函数导数(适用于Newton 迭代法) function f=nonlinerequs3(x) f=[1-0.7*cos(x(1)),0.2*sin(x(2));0.7*sin(x(1)),1+0.2*cos(x(2))]; 导数为???? ? ? ??????????????2212211 1x f x f x f x f ，对多方程则类似。 命令： fsolve(@nonlinerequs2,[0.5,0.5]) [x,n]=unmovepoints(@nonlinerequs1,[0,0],1e-6) [x,n]=newtons(@nonlinerequs2,@nonlinerequs3,[0,0],1e-6) 计算结果：（eps=0.000001） 迭代方法 X 迭代次数n 解析解 [0.52652262191818 0.50791971903685] - fsolve [0.52652266171295 0.50791973020932] - 不动点法 [0.52652130091388 0.50792028463452] 30 Newton 迭代法 [0.52652279369020 0.50791961189450] 16
C++实现 牛顿迭代 解非线性方程组
C++实现牛顿迭代解非线性方程组（二元二次为例） 求解{0=x*x-2*x-y+0.5; 0=x*x+4*y*y-4; }的方程 #include 毕业论文开题报告 信息与计算科学 非线性方程组的迭代解法 一、选题的背景和意义 =的系数矩阵具有两非线性问题是近代数学研究的主流之一，随着计算问题的日益复杂化Ax b 个明显的特点：大型化和稀疏化。大型化指系数矩阵阶数可达上万甚至更高，稀疏性指A的零元素占绝大多数对这样的A作直接三角分解，稀疏性会遭到破坏，零元素被大量填入变为非零元素，因此迫切需要新的数值方法，适用于大型稀疏线性方程，以节省储存空间和计算时间，即提高计算效 =是数值计算的重要任务，但是率，迭代法在这样的背景下得到关注和发展，求解线性方程组Ax b 大多数科学和实际问题本质上是非线性的，能做线性化的毕竟有限，对这些非线性问题是各种解决方案，常常归纳为求解一个非线性方程组，而与线性方程相比非线性方程组的求解要困难和复杂的多，计算量也大的多，现有的理论研究还比较薄弱。而对于非线性方程，一般都用迭代法求解。二、国内外研究现状、发展动态 近年来，国内外专家学者非线性方程组的迭代解法的研究兴趣与日俱增，他们多方面、多途径地对非线性方程组进行了广泛的领域性拓展（科学、物理、生产、农业等），取得了一系列研究成果。这些研究，既丰富了非线性方程组的内容，又进一步完善了非线性方程组的研究体系，同时也给出了一些新的研究方法，促进了数值计算教学研究工作的开展，推动了课程教学改革的深入进行。三、研究的主要内容，拟解决的主要问题（阐述的主要观点） 非线性的迭代法是解非线性方程组的基本途径，是数值计算中非线性方程组求根的重要工具，也是研究非线性方程组整体性质和具体分布的重要工具。就因为这样，很多专家学者对非线性方程组的迭代法进行研究。在前人研究的基础上，本文首先介绍非线性方程组迭代法的产生背景以及国内外状况，然后从数值计算的定义及理论定理出发来研究非线性方程组的迭代法的一些相关的结论，包括非线性方程组的基于不动点原理的迭代法、newton迭代法及其收敛性、非线性方程组的迭代法及其收敛性、最小二乘法、迭代法的收敛加速性等，进一步讨论非线性方程组迭代解法的收敛性质以及其他一些相关定理，以便我们更好、更清楚的看到非线性方程组和迭代法之间的联系，以及收敛和加速。 四、研究（工作）步骤、方法及措施（思路） （一）研究方法 实验五 解线性方程组的迭代法 【实验内容】 对1、设线性方程组 ?? ? ? ?? ? ? ?? ? ? ?? ? ? ??-=???????????????? ?????????????????? ? ?--------------------------211938134632312513682438100412029137264 2212341791110161035243120 536217758683233761624491131512 013012312240010563568 0000121324 10987654321x x x x x x x x x x ()T x 2,1,1,3,0,2,1,0,1,1*--= 2、设对称正定系数阵线性方程组 ?? ? ????? ??? ? ? ??---=????????????? ??????????????? ??---------------------4515229 23206019243360021411035204111443343104221812334161 2065381141402312122 00240424 87654321x x x x x x x x ()T x 2,0,1,1,2,0,1,1*--= 3、三对角形线性方程组 ?? ? ?? ? ????? ??? ? ? ??----=???????????????? ?????????????????? ??------------------5541412621357410000000014100000000141000000001410000000014100000000141000000001410000000014100000000 14100000000 1410987654321x x x x x x x x x x ()T x 1,1,0,3,2,1,0,3,1,2*---= 试分别选用Jacobi 迭代法，Gauss-Seidol 迭代法和SOR 方法计算其解。 【实验方法或步骤】 1、体会迭代法求解线性方程组，并能与消去法加以比较； 2、分别对不同精度要求，如54310,10,10---=ε由迭代次数体会该迭代法的收敛快慢； 3、对方程组2，3使用SOR 方法时，选取松弛因子ω=0.8，0.9，1，1.1，1.2等，试看对算法收敛性的影响，并能找出你所选用的松弛因子的最佳者； 4、给出各种算法的设计程序和计算结果。 程序： 用雅可比方法求的程序： function [x,n]=jacobi(A,b,x0,eps,varargin) if nargin==3 eps=1.0e-6; M=200; 线性方程组的迭代法及程序实现 学校代码:11517 学号:200810111217 HENAN INSTITUTE OF ENGINEERING 毕业论文 题目线性方程组的迭代法及程序实现 学生姓名 专业班级 学号 系 (部)数理科学系 指导教师职称 完成时间 2012年5月20日河南工程学院 毕业设计(论文)任务书 题目:线性方程组的迭代法及程序实现专业:信息与计算科学学号 : 姓名一、主要内容: 通过本课题的研究,学会如何运用有限元方法来解决线性代数方程组问题,特别是Gaussie-Seidel迭代法和Jacobi迭代法来求解线性方程组。进一步学会迭代方法的数学思想,并对程序代码进行解析与改进,这对于我们以后学习和研究实际问题具有重要的意义。本课题运用所学的数学专业知识来研究,有助于我们进一步掌握大学数学方面的知识,特别是迭代方法。通过这个课题的研究,我进一步掌握了迭代方法的思想,以及程序的解析与改进,对于今后类似实际问题的解决具有重要的意义。 二、基本要求: 学会编写规范论文,独立自主完成。 运用所学知识发现问题并分析、解决。 3.通过对相关资料的收集、整理,最终形成一篇具有自己观点的学术论文,以期能对线性方程组迭代法的研究发展有一定的实践指导意义。 4.在毕业论文工作中强化英语、计算机应用能力。 完成期限: 2012年月指导教师签名:专业负责人签名: 年月日 目录 中文摘要....................................................................................Ⅰ英文摘要 (Ⅱ) 1 综述 1 2 经典迭代法概述 3 2.1 Jacobi迭代法 3 2.2 Gauss?Seidel迭代法 4 2.3 SOR(successive over relaxation)迭代法 4 2.4 SSOR迭代法 5 2.5 收敛性分析5 2. 6 数值试验 6 3 matlab实现的两个例题8 3.1 例1 迭代法的收敛速度8 3.2 例 2 SOR迭代法松弛因子的选取 12致谢16参考文献17附录19 作业六：分别编写用Jacobi迭代法和Gauss-Seidel迭代法求解线性方程组Ax=B的标准程序，并求下列方程组的解。 可取初始向量 X(0) =(0，0，0)’； 迭代终止条件||x(k+1)-x(k)||<=10e-6 （1） = （2） = Jacobi迭代法： 流程图 开 始 判断b中的最大值 有没有比误差大 给x赋初值 进行迭代 求出x，弱到100次还没到，警告不收 结束 程序 clear;clc; A=[8,-1,1;2,10,01;1,1,-5]; b=[1;4;3]; e=1e-6; x0=[0;0;0]'; n=length(A); x=zeros(n,1); k=0; r=max(abs(b)); while r>e for i=1:n d=A(i,i); if abs(d) 程序结果（1） （2） Gauss-Seidel迭代法： 程序 clear;clc; %A=[8,-1,1;2,10,01;1,1,-5]; %b=[1;4;3]; A=[5,2,1;-1,4,2;2,-3,10]; b=[-12;20;3]; m=size(A); if m(1)~=m(2) error('矩阵A不是方阵'); end n=length(b); %初始化 N=0;%迭代次数 L=zeros(n);%分解A=D+L+U,D是对角阵，L是下三角阵，U是上三角阵U=zeros(n); D=zeros(n); G=zeros(n);%G=-inv(D+L)*U d=zeros(n,1);%d=inv(D+L)*b x=zeros(n,1); for i=1:n%初始化L和U for j=1:n if i 牛顿迭代法求解非线性方程组 非线性方程组如下： 221122121210801080 x x x x x x x ?-++=??+-+=?? 给定初值()00.0T x =，要求求解精度达到0.00001 1.首先建立函数()F X ,方程编程如下，将F.m 保存到工作路径中： function f=F(x) f(1)=x(1)^2-10*x(1)+x(2)^2+8; f(2)=x(1)*x(2)^2+x(1)-10*x(2)+8; f=[f(1),f(2)] ; 2.建立函数()DF X ,用于求方程的jacobi 矩阵,将DF.m 保存到工作路径中： function df=DF(x) df=[2*x(1)-10,2*x(2);x(2)^2+1,2*x(1)*x(2)-10]; %jacobi 矩阵是一阶偏导数以一定方式排列成的矩阵。 3.编程牛顿迭代法解非线性方程组，将newton.m 保存在工作路径中： clear,clc; x=[0,0]'; f=F(x); df=DF(x); fprintf('%d %.7f %.7f\n',0,x(1),x(2)); N=4; for i=1:N y=df\f'; x=x-y; f=F(x); df=DF(x); fprintf('%d %.7f %.7f\n',i,x(1),x(2)); if norm(y)<0.0000001 break; else end end ezplot('x^2-10*x+y^2+8',[-6,6,-6,6]); hold on ezplot('x*y^2+x-10*y+8',[-6,6,-6,6]); 运行结果如下： 0 0.0000000 0.0000000 1 0.8000000 0.8800000 2 0.9917872 0.9917117 3 0.9999752 0.9999685 求解线性方程组——超松弛迭代法 #include cout<<"输入松弛因子w (1 一、一般迭代法求解非线性方程组。 function [k,piancha,xdpiancha,xk]=diedai1(x0,k) x(1)=x0; for i=1:k x(i+1)=fun1(x(i)); piancha=abs(x(i+1)-x(i)); xdpiancha=piancha/(abs(x(i+1))+eps); i=i+1;xk=x(i); [(i-1) piancha xdpiancha xk] end if (piancha>1)&(xdpiancha>0.5)&(k>3) disp('此迭代序列发散，请重新输入新的迭代公式') return; end if (piancha<0.001)&(xdpiancha<0.0000005)&(k>3) disp('此迭代序列收敛，且收敛速度较快') return; end p=[(i-1) piancha xdpiancha xk]'; 1、function y=fun1(x) y=(10-x.^2)./2 >> [k,piancha,xdpiancha,xk]=diedai1(5,10) 此迭代序列发散，请重新输入新的迭代公式 k = 10 piancha = 2.4484e+271 xdpiancha = 1 xk = -2.4484e+271 2、function y=fun1(x) y=10./(x+2) >> [k,piancha,xdpiancha,xk]=diedai1(5,25) 此迭代序列收敛，且收敛速度较快 k = 25 piancha = 9.5676e-007 xdpiancha = 4.1300e-007 xk = 2.3166 二、第二种迭代法。 Gauss-Seidel迭代法求解线性方程组 一． 问题描述 用Gauss-Seidel 迭代法求解线性方程组 由Jacobi 迭代法中，每一次的迭代只用到前一次的迭代值。使用了两倍的存储空间，浪费了存储空间。若每一次迭代充分利用当前最新的迭代值，即在计算第i 个分量 ) 1(+k i x 时，用最新分量 ) 1(1 +k x ， ???+) 1(2 k x ) 1(1 -+k i x 代替旧分量 ) (1 k x ， ???) (2 k x ) (1 -k i x ，可以起 到节省存储空间的作用。这样就得到所谓解方程组的Gauss-Seidel 迭代法。 二． 算法设计 将A 分解成U D L A --=，则b x =A 等价于b x =--U)D (L 则Gauss-Seidel 迭代过程 ) ()1()1(k k k Ux Lx b Dx ++=++ 故 ) ()1()(k k Ux b x L D +=-+ 若设1 )(--L D 存在，则 b L D Ux L D x k k 1)(1)1()()(--+-+-= 令 b L D f U L D G 11)()(---=-=， 则Gauss-Seidel 迭代公式的矩阵形式为 f Gx x k k +=+) () 1( 其迭代格式为 T n x x x x ) ()0()0(2)0(1)0(，，，???= (初始向量)， ) (1 1 1 1 1 )()1()1(∑∑-=-+=++--=i j i i j k j ij k j ij i ii i i x a x a b a x )210i 210(n k ???=???=，，，；，，， 或者 ?? ???--=???=???==?+=∑∑-=-+=+++) (1)210i 210(111 1)()1()1()()1(i j i i j k j ij k j ij i ii i i i k i k i x a x a b a x n k k x x x ，，，；，，， 三． 程序框图 Newton迭代法求解非 线性方程 一、 Newton 迭代法概述 构造迭代函数的一条重要途径是用近似方程来代替原方程去求根。因此，如果能将非线性方程f (x )=0用线性方程去代替，那么，求近似根问题就很容易解决，而且十分方便。牛顿(Newton)法就是一种将非线性方程线化的一种方法。 设k x 是方程f (x )=0的一个近似根，把如果)(x f 在k x 处作一阶Taylor 展开，即： )x x )(x ('f )x (f )x (f k k k -+≈ (1-1) 于是我们得到如下近似方程： 0)x x )(x ('f )x (f k k k =-+ (1-2) 设0)('≠k x f ，则方程的解为： x ?=x k +f (x k ) f (x k )? (1-3) 取x ~作为原方程的新近似根1+k x ，即令： ) x ('f ) x (f x x k k k 1k -=+, k=0,1,2,… (1-4) 上式称为牛顿迭代格式。用牛顿迭代格式求方程的根的方法就称为牛顿迭代法，简称牛顿法。 牛顿法具有明显的几何意义。方程： )x x )(x ('f )x (f y k k k -+= (1-5) 是曲线)x (f y =上点))x (f ,x (k k 处的切线方程。迭代格式(1-4)就是用切线式(1-5)的零点来代替曲线的零点。正因为如此，牛顿法也称为切线法。 牛顿迭代法对单根至少是二阶局部收敛的，而对于重根是一阶局部收敛的。一般来说，牛顿法对初值0x 的要求较高，初值足够靠近*x 时才能保证收敛。若 要保证初值在较大范围内收敛，则需对)x (f 加一些条件。如果所加的条件不满足，而导致牛顿法不收敛时，则需对牛顿法作一些改时，即可以采用下面的迭代格式： ) x ('f ) x (f x x k k k 1k λ -=+, ?=,2,1,0k (1-6) 上式中，10<λ<，称为下山因子。因此，用这种方法求方程的根，也称为牛顿下山法。 牛顿法对单根收敛速度快，但每迭代一次，除需计算)x (f k 之外，还要计算 )x ('f k 的值。如果)x (f 比较复杂，计算)x ('f k 的工作量就可能比较大。为了避免计算导数值，我们可用差商来代替导数。通常用如下几种方法： 1. 割线法 如果用 1 k k 1k k x x ) x (f )x (f ----代替)x ('f k ，则得到割线法的迭代格式为： )x (f ) x (f )x (f x x x x k 1k k 1 k k k 1k --+---= (1-7) 2. 拟牛顿法 如果用 ) x (f )) x (f x (f )x (f k 1k k k ---代替)x ('f k ，则得到拟牛顿法的迭代格式为： )) x (f x (f )x (f ) x (f x x 1k k k k 2k 1k -+--- = (1-8) 3. Steffenson 法 如果用 ) x (f ) x (f ))x (f x (f k k k k -+代替)x ('f k ，则得到拟牛顿法的迭代格式为： ) x (f ))x (f x (f ) x (f x x k k k k 2k 1 k -+- =+ 第4章解线性方程组的迭代法 用迭代法求解线性方程组与第4章非线性方程求根的方法相似，对方程组进行等价变换，构造同解方程组(对可构造各种等价方程组， 如分解，可逆，则由得到)，以此构造迭代关系式 （4.1） 任取初始向量，代入迭代式中，经计算得到迭代序列。 若迭代序列收敛，设的极限为，对迭代式两边取极限 即是方程组的解，此时称迭代法收敛，否则称迭代法发散。我们将看到，不同于非线性方程的迭代方法，解线性方程组的迭代收敛与否完全决定于迭代矩阵的性质，与迭代初始值的选取无关。迭代法的优点是占有存储空间少，程序实现简单，尤其适用于大型稀疏矩阵；不尽人意之处是要面对判断迭代是否收敛和收敛速度的问题。 可以证明迭代矩阵的与谱半径是迭代收敛的充分必要条件，其中是矩阵的特征根。事实上，若为方程组的解，则有 再由迭代式可得到 由线性代数定理，的充分必要条件。 因此对迭代法(4.1)的收敛性有以下两个定理成立。 定理4.1迭代法收敛的充要条件是。 定理4.2迭代法收敛的充要条件是迭代矩阵的谱半径 因此，称谱半径小于1的矩阵为收敛矩阵。计算矩阵的谱半径，需要求解矩阵的特征值才能得到，通常这是较为繁重的工作。但是可以通过计算矩阵的范数等方法简化判断收敛的 工作。前面已经提到过，若||A||p矩阵的范数，则总有。因此，若，则必为收敛矩阵。计算矩阵的1范数和范数的方法比较简单，其中 于是，只要迭代矩阵满足或，就可以判断迭代序列 是收敛的。 要注意的是，当或时，可以有，因此不能判断迭代序列发散。 在计算中当相邻两次的向量误差的某种范数小于给定精度时，则停止迭代计算，视为方程组的近似解（有关范数的详细定义请看3.3节。） 4.1雅可比（Jacobi）迭代法 4.1.1 雅可比迭代格式 雅可比迭代计算 元线性方程组 （4.2） 写成矩阵形式为。若将式（4.2）中每个方程的留在方程左边，其余各项移到方程右边；方程两边除以则得到下列同解方程组： 记，构造迭代形式 迭代法解线性方程组作业 沈欢00986096 北京大学工学院，北京100871 2011年10月12日 摘要 由所给矩阵生成系数矩阵A和右端项b，分析系数矩阵A，并用Jacobi迭代法、GS迭代法、SOR(逐步松弛迭代法)解方程组Ax=b 1生成系数矩阵A、右端项b,并分析矩阵A 由文件”gr900900c rg.mm”得到了以.mm格式描述的系数矩阵A。A矩阵是900?900的大型稀 疏对称矩阵。于是，在matlaB中，使用”A=zeros(900,900)”语句生成900?900的零矩阵。再 按照.mm文件中的描述，分别对第i行、第j列的元素赋对应的值，就生成了系数矩阵A，并 将A存为.mat文件以便之后应用。 由于右端项是全为1的列向量，所以由语句”b=ones(900,1)”生成。 得到了矩阵A后，求其行列式，使用函数”det(A)”,求得结果为”Inf”,证明行列式太大，matlaB无法显示。由此证明，矩阵A可逆，线性方程组 Ax=b 有唯一解。 接着，判断A矩阵是否是对称矩阵(其实，这步是没有必要的，因为A矩阵本身是对称矩阵，是.mm格式中的矩阵按对称阵生成的)。如果A是对称矩阵，那么 A?A T=0 。于是，令B=A?A T,并对B求∞范数。结果显示： B ∞=0,所以，B是零矩阵，也就是：A是对称矩阵。 然后，求A的三个条件数： Cond(A)= A ? A?1 所求结果是，对应于1范数的条件数为：377.2334;对应于2范数的条件数为：194.5739;对应 于3范数的条件数为：377.2334; 1 从以上结果我们看出，A是可逆矩阵，但是A的条件数很大，所以，Ax=b有唯一解并且矩阵A相对不稳定。所以，我们可以用迭代方法来求解该线性方程组，但是由于A的条件数太大迭代次数一般而言会比较多。 2Jacobi迭代法 Jacobi迭代方法的程序流程图如图所示： 图1:Jacobi迭代方法程序流程图 在上述流程中，取x0=[1,1,...,1]T将精度设为accuracy=10?3，需要误差满足： error= x k+1?x k x k+1 数值分析实验五 班级： 10信计二班 学号：59 姓名：王志桃 分数： 一．实验名称 高斯-赛德尔迭代法解线性方程组 二．实验目的 1. 学会利用高斯赛德尔方法解线性方程组 2. 明白迭代法的原理 3. 对于大型稀疏矩阵方程组适用于迭代法比较简单 三．实验内容 利用Gauss-Seidel 迭代法求解下列方程组 ?????=++=-+=+-36123633111420238321 321321x x x x x x x x x ， 其中取→=0)0(x 。 四、算法描述 由Jacobi 迭代法中，每一次的迭代只用到前一次的迭代值，若每一次迭代充分利用当前最新的迭代值，即在计算第i 个分量)1(+k i x 时，用最新分量)1(1+k x ，???+)1(2k x )1(1-+k i x 代替旧分量)(1k x ，???)(2k x )(1-k i x ，就得到所谓解方程组的Gauss-Seidel 迭代法。 其迭代格式为 T n x x x x )()0()0(2)0(1)0(，，，???= (初始向量)， )(11111)()1( ) 1(∑∑-=-+=++--=i j i i j k j ij k j ij i ii i i x a x a b a x )210i 210(n k ???=???=，，，；，，， 或者写为 ?? ???--=???=???==?+=∑∑-=-+=+++)(1)210i 210(1111)( )1()1()()1(i j i i j k j ij k j ij i ii i i i k i k i x a x a b a x n k k x x x ，，，；，，， 五、 编码 #include 数值分析实验报告 0 a 12 K a 1,n 1 K a 2,n 1 U O M 则有： 第一步： Jacobi 迭代法 a 1n a 2n M , 则有： A D L U a n 1,n Ax b A A x D b L U (D L U)x b Dx (L U)x b x D (L U)x D b 令 J D (L U) 则称 J 为雅克比迭代矩阵 f D b 由此可得雅克比迭代的迭代格式如下： x (0) , 初始向量 x (k 1) Jx (k) f ,k 0,1,2,L 第二步 Gauss-Seidel 迭代法 Ax b (D L U )x b (D L)x Ux b x (D L) Ux (D L) b A D L U a 11 a 12 L a 1n a 11 A a 21 a 22 L a 2n a 22 M MM MO a n1 a n2 L a nn a 11 得到 D a 22 O a nn 由 a 21 0 M M O a n 1,1 a n 1,2 L 0 a nn a n1 a n2 L a n,n a 21 L M M O a n 1,1 a n 1,2 L a n1 a n2 L a n,n 1 a 12 K a 1,n 1 a 1n 0 K a 2,n 1 a 2n O M M a n 1,n 10 令 G (D L) U ，则称G 为Gauss-Seidel 迭代矩阵 f (D L) b 由此可得 Gauss-Seidel 迭代的迭代格式如下： x (0) , 初始向量 第三步 SOR 迭代法 w0 AD L U 1 ( D 1 wL ((1 w)D wU )) (D 1 wL) ((1 w)D wU ) w w w 令M w 1 (D wL), N 1 ((1 w)D wU )则有：A MN w w Ax b AM L W N M (M N )x b Mx Nx b x M Nx M b N M, 令W f Mb 带入 N 的值可有 L W ((1 w)D wU) (D wL) 1((1 w)D wU) (D wL) f 1 b w 1(D wL) 1b 1 (D wL) w 称 L W 为 SOR 迭代矩阵，由此可得 SOR 迭代的迭代格式如下： x (0) ,初始向量 二、算法程序 Jacobi 迭代法的 M 文件： function [y,n]=Jacobi(A,b,x0,eps) %************************************************* %函数名称 Jacobi 雅克比迭代函数 %参数解释 A 系数矩阵 % b 常数项 % x0 估计解向量 x (k 1) Gx (k) f ,k 0,1,2,L (k 1) f,k 0,1,2,L 华中师范大学 课程结业论文 题目：非线性方程组牛顿法及MATLAB程序 院系：数学与统计学学院 专业：数学与应用数学 年级：2014级 课堂名称：数值分析（1）实验 学生姓名：杨帅 学号：2014212643 2016年6月18 非线性方程组牛顿法及其MATLAB 程序 〔摘要〕学了《数值分析》这门课，了解到非线性方程的数值解法有：对分区间法、简单迭代法、Aitken-Steffensen 加速法、Newton 迭代法、正割法等，自然就会想到非线性方程组的数值解法有哪些呢？和非线性方程的数值解法有哪些不不同呢？ 在研究非线性方程组的数值解法之前，首先要给非线性方程组下一个合理定义；n 个变量n 个方程（n>1）的方程组表示为0 ),...,,(2 1 =n i x x x f （其中i=1,2...，n ），若i f 中至少有一个是非线性函数，则称上述 的表示为非线性方程组。在R 中记，T n x x x f f ),...,,(2 1 =，其中记 ),...,()(1n i i i x x f x f f ==且D x ∈。 若存在尣∈D ，使?(尣)=0，则称尣为非线性方程组的解。上述方程组可能有一个解或多个解，也可能有无穷多解或无解。对非线性方程组解的存在性的研究远不如线性方程组那样成熟，现有的解法也不象线性方程组那样有效。除极特殊的方程外，一般不能用直接方法求得精确解，目前主要采用迭代法求近似解。根据不同思想构造收敛于解尣的迭代序列{尣}(k=0，1，…)，即可得到求解非线性方程组的各种迭代法；但研究数学问题的时候，一般是由简单到复杂，由特殊到一般。因此要在研究非线性方程组牛顿解法的时候，首先要探究非线性方程的牛顿解法。 1.1求解线性方程组的牛顿法及其MATAB 程序 1.1.1程序设计思路 输入的量：初始值0 x 、近似根k x 的误差限tol ，近似根k x 的函数 值)(k x f 得误差限ftol ，迭代次数的最大值gxmax 、函数fnq （x ）=) (x f 数值分析方法中方程求解的直接法和迭代法 第3章 解线性方程组的直接法 一、 消元法 1. 高斯消元法(加减消元)：首先将A 化为上三角阵，再回代求解。 11121121222212n n n n nn n a a a b a a a b a a a b ?? ? ? ? ??? (1)(1)(1)(1)(1)11 121311(2)(2)(2)(2)222322 (3)(3)(3)3333()()000 00 n n n n n nn n a a a a b a a a b a a b a b ?? ? ? ? ? ? ?? ? 步骤如下： 第一步：1 11 1,2,,i a i i n a -? +=第行第行 11121121222212 n n n n nn n a a a b a a a b a a a b ?? ? ? ? ??? 111211(2)(2)(2)2222 (2)(2)(2)2 00n n n nn n a a a b a a b a a b ?? ? ? ? ??? 第二步：(2)2 (2)222,3, ,i a i i n a -?+=第行第行 111211(2)(2)(2)2222 (2)(2)(2)200n n n nn n a a a b a a b a a b ?? ? ? ? ?? ? 111213 11 (2)(2)(2)(2) 222322 (3)(3)(3) 33 33(3)(3)(3) 3 00000n n n n nn n a a a a b a a a b a a b a a b ?? ? ? ? ? ? ?? ? 类似的做下去，我们有： 第k 步：() ()k ,1, ,k ik k kk a i i k n a -?+=+第行第行。 n －1步以后，我们可以得到变换后的矩阵为：非线性方程组的迭代解法【开题报告】
常微分方程的解线性方程组的迭代法
线性方程组的迭代法及程序实现
数值分析5-用Jacobi迭代法和Gauss-Seidel迭代法求解线性方程组
牛顿迭代法求解非线性方程组的代码
求解线性方程组——超松弛迭代法(c)
十、解非线性方程(组)的迭代法和加速法
Gauss-Seidel迭代法求解线性方程组
Newton迭代法求解非线性方程
数值计算_第4章 解线性方程组的迭代法
迭代法解线性方程组
高斯-赛德尔迭代法解线性方程组精选.
实验解线性方程组的基本迭代法实验
非线性方程组牛顿迭代法(1)
解线性方程组的直接法和迭代法
- 第四章线性方程组和非线性方程组的迭代法
- 不动点迭代法和牛顿法非线性方程组求解
- 第四章线性方程组和非线性方程组的迭代法
- 解线性方程组的迭代法
- 实验八、解线性方程组的迭代法
- 线性方程组的直接法和迭代法
- 第七章解线性方程组的迭代法
- 解线性方程组的迭代法
- 不动点迭代法求解非线性方程组
- 计算方法 第八章 解线性方程组的迭代法 高斯迭代法 迭代法的收敛性
- 非线性方程组的迭代解法【开题报告】
- 第6章解线性方程组的迭代法
- 线性方程组的迭代解法
- 迭代法求解线性方程组的研究(精选.)
- 非线性方程组迭代法
- 浙大计算方法 6 迭代法的收敛性及非线性方程组的解法(6.9-10)
- 线性方程组的迭代法
- 计算方法第六章解线性方程组的迭代法
- 线性方程组的迭代法
- 解线性方程组的迭代法