1、目錄第一章 線性方程組的直接解法21.1 高斯消去法解線性方程組2 1.2 列主元高斯消去法解線性方程組3 1.3 用LU方法解線性方程組5 1.4 利用追趕法求解三對角線性方程組7第二章 解線性方程組的迭代法92.1 求解線性方程組的Jacobi迭代法9 2.2 求解線性方程組的Gauss-Seidel迭代法10 2.3 求解線性方程組的SOR迭代法12第三章 解非線性方程（組）的迭代法143.1 求解非線性方程的一般迭代法14 3.2 求解非線性方程的二分法15 3.3 求解非線性方程的Newton法17 3.4 求解非線性方程組的牛頓法19第四章 多項式插值224.1 Lagrange插

2、值多項式的構(gòu)造22 4.2 Newton插值多項式的構(gòu)造23第五章 數(shù)值積分275.1 Simpson公式和梯形公式27 5.2 復(fù)化Simpson公式和復(fù)化梯形公式27第六章 常微分方程初值問題的數(shù)值解法296.1 Euler方法和改進的Euler方法29 6.2 四階龍格庫塔方法27 第一章 線性方程組的直接解法室驗內(nèi)容簡介：如何利用幾種常見的直接解法求解線性方程組實驗?zāi)康模杭訌娋幊棠芰途幊碳记桑毩?xí)從數(shù)值分析的角度看問題。主要儀器及配套軟件：Matlab、C、Mathematica等1.1高斯消去法解線性方程組實驗項目名稱：用高斯消去法求解下列方程組Ax=b實驗內(nèi)容：題目： 已知線性方

3、程組如下：A = 2 2 2 3 2 4 1 3 9b = 1.0000 0.5000 2.5000原理：將增廣矩陣化為下三角陣，從而求解線性方程組。基本思想： 用Gauss消去法解線性方程組的基本思想是設(shè)法消去方程組的系數(shù)矩陣A的主對角線下的元素，而將線性方程組化為等價的上三角形方程組，然后再通過回代過程便可獲得原方程組的解。Matlab源代碼：function x =Gauss_solve(A, b)% x = Gauss_solve(A, b) returns the solution to the equation Ax = b,% where A is an n-by-n matri

4、x and b is a column vector of% length n (or a matrix with several such columns).n, n = size(A);for k = 1 : n-1 if A(k,k)=0 disp('A is singular,stop.') end % zero out entries of A and b using pivot A(k, k) for i = k+1 : n mik = A(i, k) / A(k, k); b(i) = b(i) - mik * b(k); for j=k+1:n A(i, j)

5、= A(i, j) - mik * A(k, j); end endend% back substitutionx=zeros(n,1);x(n)=b(n)/A(n,n);for i=n-1:-1:1 sum=0; for j=i+1:n sum=sum+A(i,j)*x(j); end x(i)=(b(i)-sum)/A(i,i);end實驗結(jié)果：在命令窗口輸入相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)Gauss_solve>>A=2 2 2;3 2 4;1 3 9A = 2 2 2 3 2 4 1 3 9>>b=1;1/2;5/2b = 1.0000 0.5000 2.5000&

6、gt;>x = Gauss_solve(A, b)x = -0.5000 1.0000 01.2 列主元高斯消去法解線性方程組實驗項目名稱：用列主元高斯消去法求解下列方程組Ax=b實驗內(nèi)容：題目： A = 0.7290 0.8100 0.9000 1.0000 1.0000 1.0000 1.3310 1.2100 1.1000b = 0.6867 0.8338 1.0000原理： 在采用Gauss消去法求解線性方程組時，小主元素可能導(dǎo)致計算失敗，故在消去法中應(yīng)避免采用絕對值很小的主元素。故在Gauss消去法中的每一步需要采取選主元的技巧。設(shè)計思想：在高斯消去法的基礎(chǔ)上，列主元消去法在每

7、次選主元時，僅依次按列選取絕對值最大的元素作為主元素，且只交換相應(yīng)的兩行（不需要進行列交換），再進行消去計算。Matlab源代碼：function x = lsolve(A, b)% x = lsolve(A, b) returns the solution to the equation Ax = b,% where A is an n-by-n matrix and b is a column vector of% length n (or a matrix with several such columns).% Gaussian elimination with partial piv

8、otingn, n = size(A);for k = 1 : n-1 % find index of largest element below diagonal in column k ik = k; for i = k+1 : n if abs(A(i, k) > abs(A(ik, k) ik = i; end end if A(ik,k)=0 disp('A is singular,stop.') end % swap with row k if ik=k for j=k:n zjz=A(k,j); A(k,j)=A(ik,j); A(ik,j)=zjz; en

9、d zjz=b(k); b(k)=b(ik); b(ik)=zjz; end % zero out entries of A and b using pivot A(k, k) for i = k+1 : n mik = A(i, k) / A(k, k); b(i) = b(i) - mik * b(k); for j=k+1:n A(i, j) = A(i, j) - mik * A(k, j); end endend% back substitutionx=zeros(n,1);x(n)=b(n)/A(n,n);for i=n-1:-1:1 sum=0; for j=i+1:n sum=

10、sum+A(i,j)*x(j); end x(i)=(b(i)-sum)/A(i,i);end實驗結(jié)果：在命令窗口輸入相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)lsolve>>A=0.729 0.810 0.900;1.000 1.000 1.000;1.331 1.210 1.100A = 0.7290 0.8100 0.9000 1.0000 1.0000 1.0000 1.3310 1.2100 1.1000>>b=0.6867;0.8338;1.0000b = 0.6867 0.8338 1.0000>>x = lsolve(A, b)x = 0.2245 0.2

11、814 0.32791.3 用LU方法解線性方程組實驗項目名稱：用LU法求解下列方程組Ax=b實驗內(nèi)容：題目： A = 2 2 2 3 2 4 1 3 9b = 1.0000 0.5000 2.5000原理： 對于線性方程組的系數(shù)矩陣A，求出單位下三角矩陣L和上三角矩陣U，使得A=LU。事實上，L的主對角線左下角元素是求解該線性方程組的Gauss消去過程中使用的乘子，并按消去過程的方式排列所得。設(shè)計思想：應(yīng)用LU分解法求解線性方程組Ax=b等價于求解兩個三角形線性方程組：（1） 求解下三角線性方程組Ly=b,得y;（2） 求解上三角線性方程組Ux=y,得x。Matlab源代碼：function

12、 L,U,x =lu_solve(A, b)% x = lu_solve(A, b) returns the solution to the equation Ax = b,% where A is an n-by-n matrix and b is a column vector of% length n (or a matrix with several such columns).% LU decomposition% n, n = size(A);% compute the first row of U and the first column of Lfor i=2:n A(i,1)

13、=A(i,1)/A(1,1);endfor k=2:n % compute the kth row of U for j=k:n for p=1:k-1 A(k,j)=A(k,j)-A(k,p)*A(p,j); end end % compute the kth column of L if kn for i=k+1:n for p=1:k-1 A(i,k)=A(i,k)-A(i,p)*A(p,k); end A(i,k)=A(i,k)/A(k,k); end endendL=zeros(n,n);for i=1:n L(i,i)=1;endfor i=2:n for j=1:i-1 L(i,

14、j)=A(i,j); endendU=zeros(n,n);for i=1:n for j=i:n U(i,j)=A(i,j); endend %solve the equations Ly=b;y=b;for i=2:n for j=1:i-1 y(i)=y(i)-L(i,j)*y(j); endend% solve the equations Ux=yx=y;x(n)=y(n)/U(n,n);for i=n-1:-1:1 for j=i+1:n y(i)=y(i)-U(i,j)*x(j); end x(i)=y(i)/U(i,i);end實驗結(jié)果：在命令窗口輸入相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)

15、lu_solve>>A=1 2 3;2 5 2;3 1 5A = 1 2 3 2 5 2 3 1 5>>b=14;18;20b = 14 18 20>>L,U,x =lu_solve(A, b)L = 1 0 0 2 1 0 3 -5 1U = 1 2 3 0 1 -4 0 0 -24x = 1 2 31.4 利用追趕法求解三對角線性方程組實驗項目名稱：用追趕法求解三對角線性方程組實驗內(nèi)容：題目：用追趕法解方程組原理：1. 消元過程：u1=c1/b1, y1=f1/b1ui=ci/(bi-ui-1*ai)yi=(fi-yi-1*ai)/(bi-ui-1*ai

16、)2. 回代過程：xn=ynxi=yi-ui*xi+1Matlab源代碼： 實驗結(jié)果： d =1.0842 -0.7368 0.6000 -0.2105 0.2421實驗注意：通過本實驗要了解追趕法求方程組的解的過程和內(nèi)涵，這種方法的消元過程和回代過程是關(guān)鍵的地方，通過這次實驗要進一步加強對對角占優(yōu)陣的了解。第二章 解線性方程組的迭代法實驗內(nèi)容簡介：如何利用幾種常見的迭代解法求解線性方程組實驗?zāi)康模壕毩?xí)引入迭代矩陣的形式來解決方程組，鞏固方程與矩陣相互關(guān)系的概念。主要儀器及配套軟件：Matlab、C、Mathematica等2.1 求解線性方程組的Jacobi迭代法實驗項目名稱：用Jacobi

17、迭代法求解方程組實驗內(nèi)容：題目：用雅可比迭代求方程組： AX=b 已知：A、b如下：A = 10 -1 2 0 -1 11 -1 3 2 -1 10 -1 0 3 -1 8b = 6 25 -11 15原理： 先找出下三角陣-L，再找出上對角陣-U，還有主對角陣D。 矩陣形式的迭代公式x=inv(D)*(L+U)*x+inv(D)*b設(shè)計思想： 利用迭代即可在循環(huán)中實現(xiàn)的原則來完成此方程組的求解，所用的迭代公式為：x=inv(D)*(L+U)*x+inv(D)*bMatlab源代碼：function iter,x=Jacobi(A,b,x0,eps,M)% x0 the initial poi

18、nt% M the maximum iteration % eps the error tolerancen=length(x0);x=x0;iter=0;for k=1:M for i=1:n sum=0; for j=1:n if j=i sum=sum+A(i,j)*x0(j);%x0 stands for x_k end end x(i)=(b(i)-sum)/A(i,i);% x stands for x_k+1 end if norm(x-x0)<eps iter=k; return end x0=x;enddisp('超過最大迭代次數(shù)')實驗結(jié)果：在命令窗口

19、輸入相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)Jacobi>>A=10 -1 2 0;-1 11 -1 3;2 -1 10 -1;0 3 -1 8A = 10 -1 2 0 -1 11 -1 3 2 -1 10 -1 0 3 -1 8>>b=6;25;-11;15b = 6 25 -11 15>>x0=0;0;0;0x0 = 0 0 0 0>>M=100M = 100>>eps=1e-4eps = 1.0000e-004>>iter,x=Jacobi(A,b,x0,eps,M)iter = 13x = 1.0000 2.0000 -1.

20、0000 1.00002. 求解線性方程組的Gauss-Seidel迭代法實驗項目名稱：用Gauss-Seidel迭代法求解方程組實驗內(nèi)容：題目：用Gauss-Seidel迭代求方程組： AX=b。方程組如2.1節(jié)中的題目。原理： 先找出下三角陣-L，再找出上對角陣-U，還有主對角陣D。 矩陣形式的迭代公式x=inv(D-L)*U*x+inv(D-L)*b設(shè)計思想： 利用迭代即可在循環(huán)中實現(xiàn)的原則來完成此方程組的求解，所用的迭代公式為：x=inv(D-L)U*x+inv(D-L)*bMatlab源代碼：function iter,x=GS(A,b,x0,eps,M)n=length(x0);i

21、ter=0;for k=1:M x=x0; % x stands for x_k, x0 stands for x_k+1 for i=1:n if i=1 sum=0; for j=2:n sum=sum+A(i,j)*x0(j); end t=(b(i)-sum)/A(i,i); x0(i)=t; end% end i=1 if 1<i<n sum=0; for j=1:i-1 sum=sum+A(i,j)*x0(j); end for j=i+1:n sum=sum+A(i,j)*x0(j); end t=(b(i)-sum)/A(i,i); x0(i)=t; end% en

22、d 1<i<n if i=n sum=0; for j=1:n-1 sum=sum+A(i,j)*x0(j); end t=(b(i)-sum)/A(i,i); x0(i)=t; end% end i=n end %end i if norm(x-x0)<eps x=x0; iter=k; return endend%end kdisp('超過最大迭代次數(shù)')實驗結(jié)果：在命令窗口輸入2.1節(jié)中一樣的相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)GS>>iter,x=GS(A,b,x0,eps,M)iter = 6x = 1.0000 2.0000 -1.0000 1

23、.00002.3 求解線性方程組的SOR迭代法實驗項目名稱：用SOR迭代法求解方程組實驗內(nèi)容：題目：用SOR迭代求方程組： AX=b。方程組如2.1節(jié)中的題目。原理： 先找出下三角陣-L，再找出上對角陣-U，還有主對角陣D。 矩陣形式的迭代公式x=inv(D-w*L)*(1-w)*D+w*U+w*inv(D-w*L)*b,其中w為松弛因子。設(shè)計思想： 利用迭代即可在循環(huán)中實現(xiàn)的原則來完成此方程組的求解，所用的迭代公式為：x=inv(D-w*L)*(1-w)*D+w*U+w*inv(D-w*L)*bMatlab源代碼：function iter,x=SOR(w,A,b,x0,eps,M)n=le

24、ngth(x0);iter=0;for k=1:M x=x0; for i=1:n if i=1 sum=0; for j=2:n sum=sum+A(i,j)*x0(j); end t=(b(i)-sum)/A(i,i); x0(i)=(1-w)*x0(i)+w*t; end% end i=1 if 1<i<n sum=0; for j=1:i-1 sum=sum+A(i,j)*x0(j); end for j=i+1:n sum=sum+A(i,j)*x0(j); end t=(b(i)-sum)/A(i,i); x0(i)=(1-w)*x0(i)+w*t; end% end

25、1<i<n if i=n sum=0; for j=1:n-1 sum=sum+A(i,j)*x0(j); end t=(b(i)-sum)/A(i,i); x0(i)=(1-w)*x0(i)+w*t; end% end i=n end %end i if norm(x0-x)<eps x=x0; iter=k; return endend%end kdisp('超過最大迭代次數(shù)')實驗結(jié)果：在命令窗口輸入2.1節(jié)中一樣的相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)SOR>>w=1.2w = 1.2000>>iter,x=SOR(w,A,b,x0,ep

26、s,M)iter = 7x = 1.0000 2.0000 -1.0000 1.0000>>w=1.5w = 1.5000>>iter,x=SOR(w,A,b,x0,eps,M)iter = 15x = 1.0000 2.0000 -1.0000 1.0000第三章 解非線性方程（組）的迭代法實驗內(nèi)容簡介：如何利用幾種常見的迭代解法求解非線性方程（組）。實驗?zāi)康模和ㄟ^本章的練習(xí)，加深理解迭代法的真諦，掌握其中的規(guī)律。主要儀器及配套軟件：Matlab、C、Mathematica等3.1 求解非線性方程的一般迭代法實驗項目名稱：通過構(gòu)造簡單迭代格式求解非線性方程實驗內(nèi)容：題

27、目：給出計算 x= 的迭代公式，討論迭代過程的收斂性，并證明 x=2。原理： 根據(jù)上式的規(guī)律可以看出，若將每個根號內(nèi)的表達式看作一個整體，則這個整體就會呈現(xiàn)出迭代的規(guī)律。設(shè)計思想：經(jīng)以上分析，令x=(2+x)0.5作為迭代公式。Matlab源代碼：disp(' 給出計算 x=(2+(2+(2+(2 + .)0.5)0.5)0.5)0.5)0.5 ')disp(' 的迭代公式，討論迭代過程的收斂性，并證明 x=2 ') fprintf('n迭代公式為 x= (2 + x)0.5 n')x1=20.5;x2=2;fazi=input('請輸入

28、精度(當(dāng)精度為0時，說明為精確值！)>>')while abs(x2-x1) > fazi x1= (2 + x1)0.5; x2= (2 + x1)0.5;endif x1 = 2 fprintf('n此時為準確值：x = %fn',x1)else fprintf('n近似值：x = %f ',x1)endfprintf('n/*/')實驗結(jié)果：3.2 求解非線性方程的二分法實驗項目名稱：利用二分法求解非線性方程實驗內(nèi)容：題目：在區(qū)間1，3用二分法求。原理：將給定的初始有根區(qū)間逐步分半，直到找到滿足一定精度要求的解。設(shè)

29、計思想：對當(dāng)前的有根區(qū)間，取其中點作為真解的近似值，若滿足精度要求則停止計算，否則，從得到的兩個小區(qū)間中選擇一個作為新的有根區(qū)間，重復(fù)上述過程。程序中重要的是循環(huán)語句的使用和有根區(qū)間端點的判斷。Matlab源代碼：function f=fun(x)f=2(-x)+exp(x)+2*cos(x)-6;function iter,x=bisection(a,b,eps)% a b are the given initial points such that f(a)*f(b)<0% eps is the error tolerancefuna=fun(a);funb=fun(b);if fu

30、na*funb>=0 disp('Error,please give the correct initial points a and b!');endk=0;% k stands for the iterationwhile 1 c=(a+b)/2; func=fun(c); k=k+1; if (abs(c-a)<=eps)|(abs(func)<=eps) x=c; iter=k; return; end if funa*func<0 b=c; else a=c; endend實驗結(jié)果：在命令窗口輸入相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)bisection&

31、gt;> a=1a = 1>> b=3b = 3>> eps=1e-4eps = 1.0000e-004>> iter,x=bisection(a,b,eps)iter = 15x = 1.8294>> fun(x)ans = 9.4896e-005下面僅給出一個基于C語言的二分法的程序供閱讀學(xué)習(xí)。#define TRUE 1/* a, c : current end points epsilon : telerance it_limit : limit of iteration number Y_a, Y_c : y values of

32、the current end points fun_f(x) : functional value at x */void main()int it_limit, it, ks; float a, b, c, epsilon, Y_a, Y_b, Y_c;float fun_f(); printf( "n Bisection Methodnn" ); while( TRUE ) printf( "Lower bound: a ? "); scanf( "%f", &a ); printf( "Upper bound

33、: c ? " ); scanf( "%f", &c ); printf( "Tolerance: epsilon ? " ); scanf( "%f", &epsilon ); printf( "Iteration limit ? " ); scanf( "%d", &it_limit ); printf( " It. a b c f(a) "); printf( " f(c) abs(f(c)-f(a)/2n" );

34、it = 0; Y_a = fun_f( a ); Y_c = fun_f( c ); if( Y_a*Y_c > 0 ) printf( " f(a)f(c) > 0 n" ); else while( TRUE ) it = it + 1; b = (a + c)/2; Y_b = fun_f( b ); printf("%3d %10.6f, %10.6f, %10.6f, %10.6f, %10.6f", it, a, b, c, Y_a, Y_c ); printf( " %12.3en", fabs(Y_c -

35、 Y_a)/2); if( it > it_limit ) break; if( fabs( b - a ) < epsilon ) break; if( Y_a*Y_b <= 0 ) c = b; Y_c = Y_b; else a = b; Y_a = Y_b; if ( it<=it_limit ) printf( " Tolerance is satisfied. n" ); if ( it>it_limit ) printf( "Iteration limit exceeded.n" ); printf( &quo

36、t;Final result: Approximate root = %g n", b ); printf("nType 1 to continue, or 0 to stop.n"); scanf("%d",&ks); if( ks != 1 ) exit(0); float fun_f(x) /* Definition of Function */float x; float f; f = x*x*x - 3*x*x - x + 3; return( f ); 3.3 求解非線性方程的Newton法實驗項目名稱：利用Newton法求

37、解非線性方程實驗內(nèi)容：題目：用牛頓法解。原理：牛頓法是一種非常重要的迭代方法，對于迭代方法而言，最重要的是迭代函數(shù)的構(gòu)造。Newton法的迭代函數(shù)的構(gòu)造依據(jù)是“由當(dāng)前點的切線來代替原曲線，將該切線與x軸的交點作為下一個迭代點”。設(shè)計思想：根據(jù)方程，給出Newton法的迭代函數(shù)，用循環(huán)語句來實現(xiàn)迭代過程。Matlab源代碼：function f=fun(x)f=exp(x)-1.5-atan(x);function f=fun_derivative(x)% give the derivative function of given functionf=exp(x)-1/(1+x2);functi

38、on iter,x=newton_equation(x0,M,eps1,eps2)% x0 the initial point% eps1 and eps2 the error tolerance% M the maximum iterationiter=0;v=fun(x0);if abs(v)<eps2 x=x0; returnendfor k=1:M x1=x0-v/fun_derivative(x0); v=fun(x1); iter=iter+1; if (abs(x1-x0)<eps1)|(abs(v)<eps2) x=x1; return end x0=x1;

39、endend實驗結(jié)果：在命令窗口輸入相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)newton_equation>>x0=-7x0 = -7>>M=20M = 20>>eps1=1e-3eps1 = 1.0000e-003>>eps2=1e-5eps2 = 1.0000e-005>>iter,x=newton_equation(x0,M,eps1,eps2)iter = 4x = -14.1011>>fun(x)ans = -7.9958e-007 下面僅給出一個基于C語言的Newton法的程序供閱讀學(xué)習(xí)。/* Newton's S

40、cheme */#include <stdio.h>#include <stdlib.h>#include <math.h>#define TRUE 1/* x : current x value (approximation for the root) xb: previous value of x y : y value for x y_derivative : derivative of y n : iteration counter */void main() int i, k, n;float epsilon, x, xb, y, y_deriva

41、tive, error;void func(); printf( "n Newton's Scheme n" ); printf( "nTolerance ? " ); scanf( "%f", &epsilon ); while( TRUE) printf( "Initial guess ? " ); scanf( "%f", &x ); error=1.0e10; xb = x; n = 0; printf( "n It.No. x(n-1) y(n-1)

42、x(n)" ); while( error>epsilon ) n = n + 1; func( x, &y, &y_derivative ); x = x - y/y_derivative; /* finds new x. */ printf( "n %3d %12.5e %12.5e %12.5e ", n, xb, y, x ); error=fabs(x-xb); xb = x; printf( "n-n" ); printf( " Final solution = %g n", x ); pr

43、intf( "-n" ); printf( "nType 1 to continue, or 0 to stop. " ); scanf( "%d", &k ); if( k != 1 ) exit(0); void func(x, y, y_derivative) /*Computes y and y_derivative */float x, *y, *y_derivative; *y = pow(x, 3) - 5.0*pow(x, 2) + 6.*x; *y_derivative = 3.0*pow(x, 2) - 1

44、0.0*x + 6.; return; 3.4 求解非線性方程組的牛頓法實驗項目名稱：利用牛頓法求解非線性方程組實驗內(nèi)容：題目：用牛頓法求解非線性方程組 該方程組可等價寫成，其中。原理： 解非線性方程組的牛頓法的迭代公式為 其中，表示的Jacobi矩陣，即設(shè)計思想：根據(jù)非線性方程組，給出牛頓法的迭代函數(shù)，用循環(huán)語句來實現(xiàn)迭代過程。Matlab源代碼：%程序funf.m給出了非線性方程組所對應(yīng)的非線性函數(shù)function f=funf(x)f=-cos(x(1)/81+x(2)2/9+sin(x(3)/3-x(1); sin(x(1)/3+cos(x(3)/3-x(2);-cos(x(1)/9+

45、x(2)/3+sin(x(3)/6-x(3);%程序發(fā)funJ.m給出了非線性函數(shù)的雅可比矩陣function J=funJ(x)J=sin(x(1)/81-1 x(2)*2/9 cos(x(3)/3; cos(x(1)/3 -1 -sin(x(3)/3;sin(x(1)/9 1/3 cos(x(3)/6-1;%程序Newton.m為解非線性方程組的牛頓法的主程序function m,r=Newton(x0,eps)% x0 為初始迭代點% eps 為誤差容限% r 最終所得的近似解% m 迭代次數(shù)n = length(x0);Fx= funf(x0);J=funJ(x0);r=x0;m=0;

46、fprintf('n 第%d次迭代所得x的值: n',m);disp(x0');for i=1:1000 d=-inv(J)*Fx; x=x0+d; x0=x; Fx=funf(x0); J=funJ(x0); fprintf('n 第%d次迭代所得x的值: n',i); disp(x0'); if norm(d)<eps r=x; m=i; return endend實驗結(jié)果：在命令窗口輸入相關(guān)數(shù)據(jù)，運行上面定義的函數(shù)Newton>> x0=0;0;0>> eps=1e-4>> m,r=Newton(x0,eps) 第0次迭代所得的x的值: 0 0 0 第1次迭代所得的x的值: -0.0129 0.3290 -0.0017 第2次迭代所得的x的值: -0.0000 0.3333 -0.0000 第3次迭代所得的x的值: 0.0000 0.3333 0.0000m = 3r = 0.0000 0