数值分析教学课件（华南理工大学）Ch10.数值算法实现.pptVIP

Ch10. 数值算法实现 §1. 线性方程组解法 1、三角形线性方程组解法 以上三角形线性方程组 为例。 回代： %文件 uptri.m function u = uptri ( a, b ) n= size(a,1) ; x(n)= b(n) / a(n, n) ; for i = n-1:-1:1 s=0 ; for j = i+1: n s=s+a(i, j) * x( j ) ; end x(i) = ( b(i) –s) / a(i, i) ; end u = x ; 2、顺序Gauss消去法 （1）消去过程： 第 k 步，计算 （2）回代过程： %文件 gauss.m function u = gauss (a, b) n = length (b) ; for k=1: n –1 for i = k+1 : n mult = a ( i, k) / a (k, k) ; for j =k +1: n % if abs ( a( k, k) ) 1e–6 a (i, j) = a (i, j) – mult * a(k, j) ; % else % disp (‘顺序Gauss消去法失败’); % pause ； % exit ； % end end b (i) = b (i) – mult * b (k) ; end end x(n)= b(n) / a(n, n) ; for i = n-1:-1:1 s=0 ; for j = i+1: n s=s+a(i, j) * x( j ) ; end x(i) = ( b(i) –s) / a(i, i) ; end u = x ; 例： % 主文件main.m a=[6, -2, 2, 4; 12, -8, 6, 10; 3,-13, 9, 3; -6, 4, 1, -18 ]; b=[16, 26, -19, -34]; x= gauss (a, b); disp ( ‘方程组解为:’ ); x 则有: main 方程组解为: x= 3 1 -2 1 3、Jacobi迭代法 % jacobi.m function y = jacobi ( a, b, x0) D = diag ( diag (a) ) ; U = - triu (a, 1) ; L = - tril (a, -1) ; B = D \ ( L+U) ; f = D \ b ; y = B* x0 + f ; n=1; while norm (y –x0 ) = 1.0e –6 x0 = y ; y = B * x0 + f ; n = n +1; end y n 例: P249 例7.21 a =[10, -1, 0; -1, 10, -2; 0, -2, 10]; b =[9; 7; 6]; jacobi ( a, b, [0; 0; 0] ) y = 0.9958 0.9579 0.7916 n = 11 §2. 方程求根 1、二分法 % erfen.m function y = erfen (fun,a, b, esp ) if nargin 4 esp =1e –4 ; end if feval (fun, a) * feval (fun, b) 0 n = 1 ; c = (a+ b) / 2 ; while (b-a) esp if feval ( fun, a) * feval (fun, c) 0 b = c ; c = ( a+b) / 2 ; elseif feval ( fun, c) * feval (fun, b) 0 a = c ; c = ( a+b) / 2 ; else y = c ; esp = 10000 ; end n= n