計(jì)算機(jī)應(yīng)用基礎(chǔ)-3-線性與非線性方程(組)求解

1、第三章第三章 代數(shù)方程（組）的代數(shù)方程（組）的Matlab求解求解3.1 線性方程（組）求解線性方程（組）求解一方程組解個(gè)數(shù)的判定一方程組解個(gè)數(shù)的判定二高斯消元法求解二高斯消元法求解三三LU矩陣分解求解矩陣分解求解四四Jacobi迭代求解迭代求解五五 調(diào)用內(nèi)部函數(shù)求解調(diào)用內(nèi)部函數(shù)求解3.2 非線性方程（組）求解非線性方程（組）求解3.1 線性方程（組）求解線性方程（組）求解一、方程組解個(gè)數(shù)及求解一、方程組解個(gè)數(shù)及求解Ax=Bmn2m1mn22221n11211aaaaaaaaaAmp2m1mp22221p11211bbbbbbbbbB對(duì)于矩陣對(duì)于矩陣3.1 線性方程（組）求解線性方程（組）求解

2、解的判定矩陣解的判定矩陣,又叫增廣矩陣又叫增廣矩陣mp2m1mmn2m1mp22221n22221p11211n11211bbbaaabbbaaabbbaaac線性方程組解的判定定理線性方程組解的判定定理: (1) 當(dāng)當(dāng) m=n 且且 rank(A)=rank(C)=n 時(shí)時(shí) 方程組有唯一解方程組有唯一解 x=A-1B x=inv(A)*B (2) 當(dāng)當(dāng)rank(A)=rank(C)=rn 時(shí)時(shí)方程組有無窮多解方程組有無窮多解:為為任任意意數(shù)數(shù)其其中中rn2 , 1i,axaxaxaxirnrn2211Z=null(A) % 求解求解A矩陣的化零矩陣矩陣的化零矩陣Z=null(A, r) %

3、求解求解A矩陣的化零矩陣的規(guī)范形式矩陣的化零矩陣的規(guī)范形式若若rank(A)norm(A*x-B)=7.47e-015 1 2 3 4 5 1 4 3 2 1 X= 4 2 1 3 2 4 3 3 4 1 3 2 2 4 A=1 2 3 4; 4 3 2 1; 1 3 2 4; 4 1 3 2; B=5 1; 4 2; 3 3; 2 4; X=inv(A)*BX1 = -9/5, 12/5 28/15, -19/15 58/15, -49/15 -32/15, 41/153.1 線性方程（組）求解線性方程（組）求解精確解：精確解： x1=inv(sym(A)*B檢驗(yàn)：檢驗(yàn)： norm(doub

4、le(A*x1-B)=0 【例例3-23-2】 A=1 2 3 4; 2 2 1 1; 2 4 6 8; 4 4 2 2 B=1; 3; 2; 6; C=A B; rank(A), rank(C)rank(A)=rank(C)=2Z=null(A,r)Z = 2.0000 3.0000 -2.5000 -3.5000 1.0000 0 0 1.0000X0= 0.9542; 0.7328; -0.0763; -0.2977X0=pinv(A)*B3.1 線性方程（組）求解線性方程（組）求解解出規(guī)范化零空間：解出規(guī)范化零空間：得出一個(gè)特解：得出一個(gè)特解：3.1 線性方程（組）求解線性方程（組）求

5、解解出全部解：解出全部解： 2 3 0.9542X = a1 -2.5 + a2 -3.5 + 0.7328 1 0 -0.0763 0 1 -0.2977驗(yàn)證得出的解：驗(yàn)證得出的解： a1=randn(1); a2=rand(1); x=a1*Z(:, 1)+a2*Z(:, 2)+x0 norm(A*x-B)【例例3-33-3】3.1 線性方程（組）求解線性方程（組）求解 1 2 3 4 1 2 2 1 1 X = 2 2 4 6 8 3 4 4 2 2 4 B=1:4; C=A B; rank(A), rank(C)最小二乘解：最小二乘解： x=pinv(A)*B檢驗(yàn)：檢驗(yàn)： norm(A

6、*x-B)=7.47e-015 二二 高斯消去法求解線性方程組高斯消去法求解線性方程組 高斯消去法它不需要方程組解的初值，也不需要重復(fù)迭代計(jì)高斯消去法它不需要方程組解的初值，也不需要重復(fù)迭代計(jì)算，它通過消去和回代兩個(gè)過程就可以直接求出方程組的算，它通過消去和回代兩個(gè)過程就可以直接求出方程組的解。和前面的迭代法一樣，在計(jì)算過程中，也有可能發(fā)散解。和前面的迭代法一樣，在計(jì)算過程中，也有可能發(fā)散而得不到方程組的解，但可以通過一些其它方法解決。而得不到方程組的解，但可以通過一些其它方法解決。設(shè)一個(gè)設(shè)一個(gè)n n元方程組，以矩陣形式表示為元方程組，以矩陣形式表示為：n21n213n2n1nn22221n1

7、1211bbb xxx aaaaaaaaa3.1 線性方程（組）求解線性方程（組）求解1 消元過程消元過程 aaaaaaaaaA3n12n11n1n212212111n11211111 n121111bbbb（1）令）令 ，進(jìn)行矩陣的代數(shù)運(yùn)算，可消去第一，進(jìn)行矩陣的代數(shù)運(yùn)算，可消去第一列的元素（列的元素（2 n) 11111 i1 iaam(2) 同理，可消去第二列的元素（同理，可消去第二列的元素（3 n)（3） 可得任意步消元公式，假設(shè)第可得任意步消元公式，假設(shè)第k步，有：步，有：3.1 線性方程（組）求解線性方程（組）求解 kkikki1kikkjikkij1kijkkkkikikbmbb

8、n, 2k, 1kj , iamaaaam（4） 經(jīng)過經(jīng)過n-1步消元后，系數(shù)矩陣稱為上三角矩陣步消元后，系數(shù)矩陣稱為上三角矩陣 a000aaaaa0aaaaAnnnn1n1nn222n , 1n1n2n222221n111n1121111n3.1 線性方程（組）求解線性方程（組）求解 nn2211nbbbb回歸過程：回歸過程：（1）首先解最后一個(gè)方程）首先解最后一個(gè)方程（2）依次回代，可得到任意公式）依次回代，可得到任意公式: 1, 3n, 2nk,a/xabxkkkn1kjjkkjkkk3.1 線性方程（組）求解線性方程（組）求解高斯消元法通用程序高斯消元法通用程序function x=U

9、pTriangleFun(A,b)% 求上三角系數(shù)矩陣的線性方程組求上三角系數(shù)矩陣的線性方程組Ax=b的解的解% A為線性方程組的系數(shù)矩陣為線性方程組的系數(shù)矩陣% b為線性方程組中的常數(shù)向量為線性方程組中的常數(shù)向量% x為線性方程組中的解為線性方程組中的解n=size(A);N=n(1);for k=N:-1:1 if(k1) s=A(k,1:(k-1)*x(1:(k-1),1); else s=0; end x(k,1)=(b(k)-s)/A(k,k);end3.1 線性方程（組）求解線性方程（組）求解高斯全主元消去法高斯全主元消去法對(duì)于求解精度較高，或方程數(shù)量較多，易采用全主元法對(duì)于求解精

10、度較高，或方程數(shù)量較多，易采用全主元法n21n213n2n1nn22221n11211bbb xxx aaaaaaaaa其增廣矩陣為其增廣矩陣為(A(1)|b(1),尋找尋找A(1)中絕對(duì)值最大為主元素中絕對(duì)值最大為主元素 0amaxa1j , inj1ni111j ,1i將此元素進(jìn)行行列交換，換為第一個(gè)元素，再進(jìn)行消將此元素進(jìn)行行列交換，換為第一個(gè)元素，再進(jìn)行消元過程元過程3.1 線性方程（組）求解線性方程（組）求解對(duì)于任意的第對(duì)于任意的第k-1次，其最大的主元素為：次，其最大的主元素為： 0amaxakj , injknikkkj ,ki進(jìn)過進(jìn)過n-1次主元消去之后，對(duì)后得到上三角矩陣。次

11、主元消去之后，對(duì)后得到上三角矩陣。通用算法通用算法3.1 線性方程（組）求解線性方程（組）求解function x,xA= AllGaussFun(A,b) % AllGaussFun.m為用高斯全主元消去法求解線性方程組的解為用高斯全主元消去法求解線性方程組的解% A為線性方程組的系數(shù)矩陣為線性方程組的系數(shù)矩陣% b為線性方程組的常數(shù)向量為線性方程組的常數(shù)向量% x為線性方程組的解為線性方程組的解% xA為消元后的系數(shù)矩陣為消元后的系數(shù)矩陣N = size(A);n = N(1);indexl = 0;indexr = 0;order = 1:n; %記錄未知數(shù)順序的向量記錄未知數(shù)順序的向量

12、for i=1:(n-1) me = max(max(abs(A(i:n,i:n); %選取全主元選取全主元 for k=i:n for r=i:n if(abs(A(k,r)=me) indexl = k; indexr = r; k=n; break; end end endtemp = A(i,1:n); A(i,1:n) = A(indexl,1:n); A(indexl,1:n) = temp; bb = b(indexl); b(indexl)=b(i); b(i) = bb; %交換主行交換主行 temp = A(1:n,i); A(1:n,i) = A(1:n,indexr);

13、 A(1:n,indexr) = temp; %交換主列交換主列 pos = order(i); order(i) = order(indexr); order(indexr) = pos; %主列的交換會(huì)造成未知數(shù)順序的變化主列的交換會(huì)造成未知數(shù)順序的變化 for j=(i+1):n if(A(i,i)=0) disp(對(duì)角元素為對(duì)角元素為0！); return; end l = A(j,i); m = A(i,i); A(j,1:n)=A(j,1:n)-l*A(i,1:n)/m; b(j)=b(j)-l*b(i)/m; endendx=UpTriangleFun(A,b);y=zeros(

14、n,1);for i=1:n for j=1:n if(order(j)=i) y(i)=x(j); end endend x=y;xA =A;例例【3-4】求多元線性方程組的解求多元線性方程組的解16x410 x3x21x3xx56xx2xx34434314321求解程序求解程序“l(fā)i3_4.m”x =（2.40, -19.20, -1.00, 4.00)3.1 線性方程（組）求解線性方程（組）求解三三 LU矩陣分解求解矩陣分解求解線性方程組變?yōu)椋壕€性方程組變?yōu)椋篖Ux=b令：令：Y=Ux則：則：LY=b3.1 線性方程（組）求解線性方程（組）求解解兩個(gè)三角行方程組，解兩個(gè)三角行方程組，可得

15、解?？傻媒?。A=LU1ll1l1L2n1n21nnn222n11211uuuuuuU其中其中function x=LUfun(A,b)% 實(shí)現(xiàn)實(shí)現(xiàn)LU分解求線性方程組分解求線性方程組Ax=b的解的解% A為線性方程組的系數(shù)矩陣為線性方程組的系數(shù)矩陣% b為線性方程組中的常數(shù)向量為線性方程組中的常數(shù)向量% x為線性方程組中的解為線性方程組中的解r=rank(A);% flag=isexist(A,b) %判斷方程組解的情況判斷方程組解的情況if flag=0 disp(該方程組無解！該方程組無解！); x=; return;else m,n=size(A); L,U,P=lu(A); %lu函數(shù)

16、為函數(shù)為MATLAB提供的函數(shù)實(shí)現(xiàn)矩陣提供的函數(shù)實(shí)現(xiàn)矩陣LU分分解解 b=P*b; % 解解Ly=b y(1)=b(1); if m1 for k=2:m y(k)=b(k)-L(k,1:k-1)*y(1:k-1); end end y=y; % 解解Ux=y得原方程組的一個(gè)特解得原方程組的一個(gè)特解 x0(r)=y(r)/U(r,r); if r1 for k=r-1:-1:1; x0(k)=(y(k)-U(k,k+1:r)*x0(k+1:r)/U(k,k); end end x0=x0; % 如果方程組有唯一解如果方程組有唯一解 if flag=1 x=x0; return; else %如

17、果方程組有無窮多解如果方程組有無窮多解 format rat; z=null(A,r); %求出對(duì)應(yīng)齊次方程組的基礎(chǔ)解系求出對(duì)應(yīng)齊次方程組的基礎(chǔ)解系 mz,nz=size(z); x0(r+1:n)=0; for k=1:nz t=sym(char(107 48+k); K(k)=t; %取取K=k1,k2,.; end x=x0; for k=1:nz x=x+K(k)*z(:,k); %將方程組的通解表示為特解加對(duì)應(yīng)齊次通解形式將方程組的通解表示為特解加對(duì)應(yīng)齊次通解形式 end endend 例例 利用利用LU分解求解如下線性方程組的解分解求解如下線性方程組的解4362x-3522532x

18、1753x43214343214321xxxxxxxxxxx3.1 線性方程（組）求解線性方程（組）求解四四 Jacobi迭代法求解迭代法求解迭代法是對(duì)任意給定的初始向量，按照某種方法逐步生迭代法是對(duì)任意給定的初始向量，按照某種方法逐步生成近似解序列，是解序列的極限為方程組的解。成近似解序列，是解序列的極限為方程組的解。設(shè)線性方程組：設(shè)線性方程組： Ax=b若若A非奇異，非奇異，b不等于不等于0，有唯一解，有唯一解x*，上式變換為：，上式變換為：x=Tx+C取一初始向量取一初始向量x(0)為上述方程的近似解為上述方程的近似解x（k+1）=Tx(k)+C滿足條件：滿足條件： *kkxxlim3.

19、1 線性方程（組）求解線性方程（組）求解此迭代法收斂，此迭代法收斂，x*為原方程組的解。為原方程組的解。Jacobi迭代法迭代法 工程常用此方法解決稀疏矩陣的迭代問題。工程常用此方法解決稀疏矩陣的迭代問題。n21n21nn2n1nn22221n11211bbb xxx aaaaaaaaa系數(shù)矩陣為系數(shù)矩陣為 aaaaaaaaaAnn2n1nn22221n112113.1 線性方程（組）求解線性方程（組）求解若若A非奇異，方程組非奇異，方程組Ax=b有唯一解有唯一解x=A-1b.方程組通過變換為：方程組通過變換為：nnn1nnn1n,n1nn1nn222n22n2122212111n11n121

20、1121atxaaxaax atxaaxaaxatxaaxaax令令iiiiiiijijatgaab3.1 線性方程（組）求解線性方程（組）求解nnnn21nn2nn222121nn12121gxbxbx gxbxbxgxbxbx aaaDnn2211令令有：有：bDg,ADIADDB1113.1 線性方程（組）求解線性方程（組）求解gBxx算法：算法： （1）設(shè)定一個(gè)初始向量）設(shè)定一個(gè)初始向量x(0), 將初始值帶入矩陣方程，得到新將初始值帶入矩陣方程，得到新向量向量x(1）。）。 （2）將新向量）將新向量x（1）帶入矩陣方程，得到另一個(gè)向量帶入矩陣方程，得到另一個(gè)向量x(2)。 (3)對(duì)于

21、任意對(duì)于任意k步迭代，迭代格式為：步迭代，迭代格式為： gBxxk1k如果序列如果序列x(k)收斂于收斂于x*，則說明，則說明x*為方程組的解為方程組的解Jacobi迭代程序：迭代程序：3.1 線性方程（組）求解線性方程（組）求解function x,n=Jacobifun(A,b,x0,eps,var)% Jacobifun為編寫的雅可比迭代函數(shù)為編寫的雅可比迭代函數(shù)% A為線性方程組的系數(shù)矩陣為線性方程組的系數(shù)矩陣% b為線性方程組的常數(shù)向量為線性方程組的常數(shù)向量% x0為迭代初始向量為迭代初始向量% eps為解的精度為解的精度% var為迭代步數(shù)控制為迭代步數(shù)控制% x為線性方程組的解為

22、線性方程組的解% n為求出所需精度的解實(shí)際的迭代步數(shù)為求出所需精度的解實(shí)際的迭代步數(shù)if nargin=3 eps= 1.0e-6; M = 200;elseif nargin=eps x0=x; x=B*x0+f; n=n+1; if(n=M) disp(Warning: 迭代次數(shù)太多，不收斂！迭代次數(shù)太多，不收斂！); return; endend3.1 線性方程（組）求解線性方程（組）求解例例【3-5】用用Jacobi迭代法計(jì)算方程組迭代法計(jì)算方程組8x5x2x2xx3x212xxx4321321321設(shè)初始值為設(shè)初始值為x0=0, 0, 0計(jì)算程序計(jì)算程序“l(fā)i3_5.m”X=2.39

23、, 2.18, -0.253.1 線性方程（組）求解線性方程（組）求解五五 調(diào)用調(diào)用Matlab內(nèi)部函數(shù)內(nèi)部函數(shù)5.1 共軛梯度法，通過最優(yōu)化求解，適用于系數(shù)矩陣為共軛梯度法，通過最優(yōu)化求解，適用于系數(shù)矩陣為方陣方陣m, 調(diào)用格式如下：調(diào)用格式如下：x, flag, relres, iter=bicg(A, b, tol)方方程程組組的的解解方方程程組組成成功功的的標(biāo)標(biāo)志志解解的的相相對(duì)對(duì)誤誤差差迭迭代代的的次次數(shù)數(shù)求求解解器器系系數(shù)數(shù)矩矩陣陣誤誤差差限限常常數(shù)數(shù)向向量量3.1 線性方程（組）求解線性方程（組）求解8x5x2x2xx3x212xxx4321321321 A=4, 1, -1;

24、2, -3, 1; 1, 2, -5; b=12, -2, 8; tol=1e-7; x, flag, relres, iter=bicg(A, b, tol)x =2.3929，2.1786 -0.2500flag = 0 relres = 1.1384e-014 iter = 33.1 線性方程（組）求解線性方程（組）求解例例3-6對(duì)于方程：對(duì)于方程：Ax+xAT=-C可采用可采用Lyapunov方法求解，其調(diào)用格式：方法求解，其調(diào)用格式： x=lyap(A,C)5.2 lyap函數(shù)函數(shù)5.3 minres 函數(shù)函數(shù) 利用殘差法求解線性方程組，其調(diào)用格式為利用殘差法求解線性方程組，其調(diào)用格

25、式為x, flag, relres, iter=minres(A, b, tol)3.1 線性方程（組）求解線性方程（組）求解3.2 非線性方程（組）求解非線性方程（組）求解一、對(duì)分法一、對(duì)分法二、直接迭代法二、直接迭代法三、牛頓迭代法三、牛頓迭代法四、割線法四、割線法五、牛頓迭代法解方程組五、牛頓迭代法解方程組六、符號(hào)求解（內(nèi)部函數(shù)）六、符號(hào)求解（內(nèi)部函數(shù)）七、非線性方程組的應(yīng)用七、非線性方程組的應(yīng)用一一 對(duì)分法對(duì)分法 對(duì)分法或稱二分法是求方程近似解的一對(duì)分法或稱二分法是求方程近似解的一種簡(jiǎn)單直觀的方法。設(shè)函數(shù)種簡(jiǎn)單直觀的方法。設(shè)函數(shù)f(x)在在a,b上連續(xù)，上連續(xù)，且且f(a)f(b)0，

26、則，則f(x)在在a,b上至少有一零點(diǎn)，上至少有一零點(diǎn)，這是微積分中的介值定理，也是使用對(duì)分法的這是微積分中的介值定理，也是使用對(duì)分法的前提條件。計(jì)算中通過對(duì)分區(qū)間，逐步縮小區(qū)前提條件。計(jì)算中通過對(duì)分區(qū)間，逐步縮小區(qū)間范圍的步驟搜索零點(diǎn)的位置。間范圍的步驟搜索零點(diǎn)的位置。 如果我們所要求解的方程從物理意義上來如果我們所要求解的方程從物理意義上來講確實(shí)存在實(shí)根，但又不滿足講確實(shí)存在實(shí)根，但又不滿足f(a)f(b)0，這，這時(shí)，我們必須通過改變時(shí)，我們必須通過改變a和和b的值來滿足二分法的值來滿足二分法的應(yīng)用條件。的應(yīng)用條件。 3.2 非線性方程（組）求解非線性方程（組）求解計(jì)算計(jì)算f(x)=0的

27、一般計(jì)算步驟如下：的一般計(jì)算步驟如下： 1、輸入求根區(qū)間、輸入求根區(qū)間a,b和誤差控制量和誤差控制量，定義函數(shù)，定義函數(shù)f(x)。 2、判斷、判斷: 如果如果f(a)f(b)0則轉(zhuǎn)下，否則，重新輸入則轉(zhuǎn)下，否則，重新輸入a和和b的值。的值。 3、計(jì)算中點(diǎn)、計(jì)算中點(diǎn) x=(a+b)/2以及以及f(x)的值的值 分情況處理分情況處理（1）|f(x)|：停止計(jì)算：停止計(jì)算x*=x，轉(zhuǎn)向步驟，轉(zhuǎn)向步驟4（2）f(a)f(x)0：修正區(qū)間：修正區(qū)間a,xa,b，重復(fù)，重復(fù)3（3）f(x)f(b)0) t=(a+b)/2;r=compute_bisect(f,t,b,tol); else if(fa*fa

28、b=0) r=(a+b)/2; else if(abs(b-a)clear all;tol=1e-6;g=bisectfun(x3-3*x+1, -1, 1,tol)g = 0.34733.2 非線性方程（組）求解非線性方程（組）求解二二 直接迭代法直接迭代法 對(duì)給定的方程對(duì)給定的方程f(x)=0，將它轉(zhuǎn)換成等價(jià)形式：，將它轉(zhuǎn)換成等價(jià)形式： 。給定初。給定初值值x0，由此來構(gòu)造迭代序列，由此來構(gòu)造迭代序列 k=1,2,，如果迭代收斂，即，如果迭代收斂，即 有有 ,則就是方程則就是方程f(x)=0的根。在計(jì)算中當(dāng)?shù)母Ｔ谟?jì)算中當(dāng) 小于給定的小于給定的精度控制量時(shí)，取精度控制量時(shí)，取 為方程的根。

29、為方程的根。 對(duì)于方程對(duì)于方程 構(gòu)造的多種迭代格式構(gòu)造的多種迭代格式 ， 怎樣判斷構(gòu)怎樣判斷構(gòu) 造的迭代格式是否收斂？收斂是否與迭代的初值有關(guān)？根據(jù)數(shù)造的迭代格式是否收斂？收斂是否與迭代的初值有關(guān)？根據(jù)數(shù)學(xué)知識(shí)，我們可以直接利用以下收斂條件：學(xué)知識(shí)，我們可以直接利用以下收斂條件： )(xx)x(xk1kbxxkkkk)(limlim1)b(bkkxx11kxb0)(xf3.2 非線性方程（組）求解非線性方程（組）求解 2 在在a,b上可導(dǎo)，并且存在正數(shù)上可導(dǎo)，并且存在正數(shù)L1，使任意，使任意的的 ，有，有 則在則在a,b上有唯一的點(diǎn)上有唯一的點(diǎn) 滿滿足足 ， 稱稱 為為 的不動(dòng)點(diǎn)。而且迭代格式

30、對(duì)任的不動(dòng)點(diǎn)。而且迭代格式對(duì)任意初值均收斂于的不動(dòng)點(diǎn)，并有下面誤差估計(jì)式：意初值均收斂于的不動(dòng)點(diǎn)，并有下面誤差估計(jì)式：,baxa)(xb)(x, baxLx | )(|1、 當(dāng)當(dāng) 有有)(*xx*x)(x|1|01*xxLLxxkk要構(gòu)造滿足收斂條件的等價(jià)形式一般比較困難。事實(shí)上，要構(gòu)造滿足收斂條件的等價(jià)形式一般比較困難。事實(shí)上，如果如果 為為f(x)的零點(diǎn)，若能構(gòu)造等價(jià)形式的零點(diǎn)，若能構(gòu)造等價(jià)形式 ，而而 ，由，由 的邊疆性，一定存在的鄰域的邊疆性，一定存在的鄰域 ，其，其上有上有 ,這時(shí)若初值這時(shí)若初值 迭代也就收斂了。由此構(gòu)造迭代也就收斂了。由此構(gòu)造收斂迭代格式有兩個(gè)要素，其一，等價(jià)形式

31、收斂迭代格式有兩個(gè)要素，其一，等價(jià)形式 應(yīng)滿足；應(yīng)滿足； 其二，初值必須取自其二，初值必須取自 的充分小鄰域，其大小決定的充分小鄰域，其大小決定于函數(shù)于函數(shù)f(x)，及做出的等價(jià)形式，及做出的等價(jià)形式)(xx*x1| )(|*x)(x,*xx,*0 xxx1| )(| Lx)(xx1| )(|*x*x)(xx3.2 非線性方程（組）求解非線性方程（組）求解上述算法的通用程序?yàn)樯鲜鏊惴ǖ耐ㄓ贸绦驗(yàn)椤癴ixptfun.m”迭代法的算法：迭代法的算法： 1. 取初始點(diǎn)取初始點(diǎn)x0，最大迭代次數(shù)，最大迭代次數(shù)N和精度要求和精度要求 ， 設(shè)設(shè)k=0; 2. 計(jì)算計(jì)算xk+1= (xk); 3. 若若 ,

32、則停止計(jì)算則停止計(jì)算 4. 若若k=N, 則停止計(jì)算，否則則停止計(jì)算，否則k=k+1,轉(zhuǎn)步驟轉(zhuǎn)步驟2k1kxx3.2 非線性方程（組）求解非線性方程（組）求解function x,k,e=fixptfun(f,x0,tol,max1)X(1)=x0;for k=2:max1 X(k)=feval(f,X(k-1); k, e=abs(X(k)-X(k-1) x=X(k); if(etol) break; end if k=max1 disp(超過最大迭代次數(shù)！超過最大迭代次數(shù)！); endendx=X;3.2 非線性方程（組）求解非線性方程（組）求解例例【3-8】利用迭代法求非線性方程組在利用

33、迭代法求非線性方程組在x0=0.4的近似的近似解，誤差為解，誤差為1e-6，最大迭代次數(shù)為，最大迭代次數(shù)為30。0 xxsin采用程序采用程序”ME_31_1.m”,計(jì)算計(jì)算X=0.8767,迭代次數(shù)為迭代次數(shù)為12.3.2 非線性方程（組）求解非線性方程（組）求解三三 牛頓迭代法牛頓迭代法對(duì)方程對(duì)方程f(x)=0可構(gòu)造多種迭代格式可構(gòu)造多種迭代格式 ，牛頓迭代法是，牛頓迭代法是借助于對(duì)函數(shù)借助于對(duì)函數(shù)f(x)=0的泰勒展開而得到的一種迭代格式。的泰勒展開而得到的一種迭代格式。 將將f(x)=0在初始值在初始值x0做泰勒展開得：做泰勒展開得：)(1kkxx)(! 2)()()()(00 000

34、 xxxfxxxfxfxf取展開式的線性部分作為的近似值，則有：取展開式的線性部分作為的近似值，則有： 0 xf)xx)(x(f)x(f000設(shè)設(shè)f(x0)K0,)()(0001xfxfxx類似地，再將類似地，再將f(x)=0在在x1作泰勒展開并取其線性部分得到作泰勒展開并取其線性部分得到：)()(1112xfxfxx3.2 非線性方程（組）求解非線性方程（組）求解得到牛頓法的迭代格式得到牛頓法的迭代格式, 2 , 1,)()(1kxfxfxxkkkk牛頓法的算法：牛頓法的算法： 1. 取初始點(diǎn)取初始點(diǎn)x0，最大迭代次數(shù)，最大迭代次數(shù)N和精度要求和精度要求 ，設(shè)，設(shè)k=0; 2. 如果如果f(

35、x0)=0，則停止計(jì)算，否則計(jì)算：，則停止計(jì)算，否則計(jì)算： 3. 若若 ,則停止計(jì)算則停止計(jì)算 4. 若若k=N, 則停止計(jì)算，否則則停止計(jì)算，否則k=k+1,轉(zhuǎn)步驟轉(zhuǎn)步驟2k1kxx)x(f)x(fxxkkk1k3.2 非線性方程（組）求解非線性方程（組）求解function x,e,k,y=newtonfun(f,df,x0,delta,max1)% newtonfun為編寫利用牛頓法求解非線性方程組的根的函數(shù)為編寫利用牛頓法求解非線性方程組的根的函數(shù)% f為非線性函數(shù)；為非線性函數(shù)； % df為為f的微商；的微商；% x0為初始值為初始值% delta為給定允許的誤差；為給定允許的誤差；

36、 % max1為迭代的最大次數(shù)為迭代的最大次數(shù)% x為用牛頓法求得的非線性方程的近似解為用牛頓法求得的非線性方程的近似解% e為為x0的誤差估計(jì)；的誤差估計(jì)； % k為迭代次數(shù)；為迭代次數(shù)；% y=f(x0)x0,feval(f,x0)for k=1:max1 x1=x0-feval(f,x0)/feval(df,x0); e=abs(x1-x0); x0=x1; x1,e,k,y=feval(f,x1) if(etol) n=n+1; x1=r; fx=subs(sym(f),x1); s=fx*fa; if(s0) r=b-(x1-b)*fb/(fx-fb); else r=a-(x1-a

37、)*fa/(fx-fa); end err=abs(r-x1); endend% secantfun.m為編寫的割線法為編寫的割線法計(jì)算非線性方程的函數(shù)計(jì)算非線性方程的函數(shù)% f為非線性方程為非線性方程% a為非線性方程的左區(qū)間為非線性方程的左區(qū)間% b為非線性方程的右區(qū)間為非線性方程的右區(qū)間% tol為誤差容限為誤差容限,當(dāng)誤差容限缺當(dāng)誤差容限缺省時(shí)省時(shí) 默認(rèn)為默認(rèn)為10的的-5次方次方% r為求得的非線性方程的近似為求得的非線性方程的近似根根 % n為迭代的次數(shù)為迭代的次數(shù)3.2 非線性方程（組）求解非線性方程（組）求解例例【3-10】利用割線法求解非線性方程在利用割線法求解非線性方程在【

38、-1， 2】上的根，精度為上的根，精度為1e-6。06x8e2x2計(jì)算程序計(jì)算程序“ME_31_3.m”，結(jié)果：，結(jié)果： r=1.6731 n=183.2 非線性方程（組）求解非線性方程（組）求解非線性方程組解法大多來自于非線性方程，只不過計(jì)算非線性方程組解法大多來自于非線性方程，只不過計(jì)算換成了矩陣，其中最重要的就是雅克比矩陣。非線性方換成了矩陣，其中最重要的就是雅克比矩陣。非線性方程組一般形式如下：程組一般形式如下：0)x,x,x(f0)x,x,x(f0)x,x,x(fn21nn212n211五五 迭代求解非線性方程組迭代求解非線性方程組為了方便，記為了方便，記Tn21Tn21)f,f,f

39、(F)x,x,x(x則，非線性方程組為則，非線性方程組為 F(x)=03.2 非線性方程（組）求解非線性方程（組）求解 非線性方程組可構(gòu)造等價(jià)形式非線性方程組可構(gòu)造等價(jià)形式)x,x,x(x)x,x,x(x)x,x,x(xn21nnn2122n2111設(shè)上述函數(shù)及其一階導(dǎo)數(shù)在設(shè)上述函數(shù)及其一階導(dǎo)數(shù)在x*區(qū)域連續(xù)，且滿足區(qū)域連續(xù)，且滿足（1）初始值）初始值x(0)足夠接近足夠接近x*;(2) 對(duì)于任意在對(duì)于任意在x*領(lǐng)域的領(lǐng)域的x點(diǎn)，點(diǎn)， 為為 F(x)的的Jacobi矩陣，即矩陣，即 1xF xF3.2 非線性方程（組）求解非線性方程（組）求解 nn1nn111xfxfxfxfxF則迭代序列則迭

40、代序列xk均收斂到均收斂到x*5.1 牛頓迭代法牛頓迭代法3.2 非線性方程（組）求解非線性方程（組）求解假定假定x*為方程組的近似解，牛頓迭代法是借助于多元函數(shù)的泰為方程組的近似解，牛頓迭代法是借助于多元函數(shù)的泰勒展開而得到的一種迭代格式。勒展開而得到的一種迭代格式。 假設(shè)它的假設(shè)它的k次近似解為次近似解為x(k)，其，其泰勒展開得：泰勒展開得： 0 xxRlim0 xxRxxxFxFxF0 xk*k*kk*其中其中當(dāng)當(dāng)x(k)接近性接近性x*，可略去，可略去 從而得到線性方程從而得到線性方程 k*xxR kkxFxxF3.2 非線性方程（組）求解非線性方程（組）求解上述方程解上述方程解 x

41、近似誤差近似誤差 k*xx 1 , 0kxFxFxxxxxxxk1kkk1kk1k有有記記實(shí)際迭代時(shí)，采用如下的方法實(shí)際迭代時(shí)，采用如下的方法 kkkk1kxFxxFxxx3.2 非線性方程（組）求解非線性方程（組）求解牛頓迭代法的算法：牛頓迭代法的算法： 1. 取初始點(diǎn)取初始點(diǎn)x0，最大迭代次數(shù)，最大迭代次數(shù)N和精度要求和精度要求 ，設(shè)，設(shè)k=0; 2. 計(jì)算計(jì)算F(x(k)和和F(x(k）)； 3. 若若Jacobi矩陣矩陣F(x(k）)奇異，則停止計(jì)算奇異，則停止計(jì)算 ,否則計(jì)算否則計(jì)算 求得求得 x（k);令令 4若若 ，停止計(jì)算，停止計(jì)算 4. 若若k=N, 則停止計(jì)算，否則則停止計(jì)

42、算，否則k=k+1,轉(zhuǎn)步驟轉(zhuǎn)步驟2k1kxx kkxFxxF k1kxxx3.2 非線性方程（組）求解非線性方程（組）求解function x,n,err=NewtonGroudfun(f,jf,x0,tol,max1)x=x0;y=feval(f,x0);for k=1:max1 j=feval(jf,x0); q=x0-(jy); y=feval(f,q); err=norm(q-x0); x0=q; x=x,x0; n=k; if(erreps fx = subs(F,findsym(F),x0); J = zeros(m,m); for i=1:m x1 = x0; x1(i) = x

43、1(i)+h; J(:,i) = (subs(F,findsym(F),x1)-fx)/h; end % FastDownGroudfun.m為編寫利用最速下降為編寫利用最速下降算法求非線性方程組的解的函數(shù)算法求非線性方程組的解的函數(shù)% F為非線性方程組為非線性方程組% x0為給定的初始值為給定的初始值% h數(shù)值微分增量步數(shù)值微分增量步% eps為解的精度為解的精度% r為求得非線性方程組的一組解為求得非線性方程組的一組解% n為迭代步數(shù)為迭代步數(shù)3.2 非線性方程（組）求解非線性方程（組）求解lamda = fx/sum(diag(transpose(J)*J); r=x0-J*lamda;

44、 fr = subs(F,findsym(F),r); tol=dot(fr,fr); x0 = r; n=n+1; if(n100000) disp(迭代步數(shù)太多，可能不收斂！迭代步數(shù)太多，可能不收斂！); return; endendformat short;3.2 非線性方程（組）求解非線性方程（組）求解六六 Matlab內(nèi)部函數(shù)（符號(hào)函數(shù)）內(nèi)部函數(shù)（符號(hào)函數(shù)）1. 內(nèi)部符號(hào)函數(shù)內(nèi)部符號(hào)函數(shù)fzero()可以求解非線性方程，其調(diào)用可以求解非線性方程，其調(diào)用格式如下：格式如下：2, 1, 0,zeroppoptionxfunfx 例例【3-11】求解非線性方程在求解非線性方程在x0=1 時(shí)解時(shí)解求解程序求解程序”ME_31_9.m” X=-3.043.2 非線性方程（組）求解非線性方程（組）求解3sinexx2. 內(nèi)部符號(hào)函數(shù)內(nèi)部符號(hào)函數(shù)fsolve()可以求解非線性方程組，采可以求解非線性方程組，采用的是最小二乘法，其調(diào)用格式如下：用的是最小二乘法，其調(diào)用格式如下：2p, 1p,option, 0 x,funfslovex 例例【3-11】求解非線性方程組在求解非線性方程組在x0=1 1時(shí)解時(shí)解0 xx3xx20 xxx4x2122221211求解程序求解程序”ME_31_6.m” X=0.25 03.2 非線性方程（組）求解非線性方程（組）求解3.2 非線性方