matlab数值分析03-方程组求解.ppt

线性方程组的计算机求解 矩阵方程组的计算机求解 矩阵方程组的计算机求解 1 线性方程组的计算机求解 Ax=B （1）若方程组有唯一解。则用Matlab语言得到该方程组的解x=inv(A)*B.但是inv（）函数的调用也有值得注意之处。例如：若A为奇异或接近奇异的则利用此函数有可能产生错误的结果。 若采用符号运算工具箱，则可以直接使用inv（）函数，如果能得出方程的解，则解是唯一的，如果出现错误信息，再考虑其它的情形。 采用如下方式可求得精确的解： x1=inv(sym(A))*B x1 = [ -9/5, 12/5] [ 28/15, -19/15] [ 58/15, -49/15] [ -32/15, 41/15] norm(double(A*x1-B)) ans = 0 （2）若方程组有无穷多解 （齐次方程） 在Matlab语言中可以由null（）直接求出，其调用格式为： Z=null（A） 求解A的化0矩阵 Z=null（A，‘r’） 求解A的化0矩阵的规范形式 注：null（）函数也可用于符号变量描述方程的解析解问题，其中Z的列数为n-r，而各列构成的向量又称为矩阵A的基础解系。 求解非齐次方程组也是比较简单的，只要能求出该方程的任意一个特解x0 则原非齐次方程组的解为x=x1+x0。其实，在Matlab中求解该方程的一个特解并非难事，用 x0=pinv（A）*B 即可求出 clear all; format long 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)] Z=null(A,r) x0=pinv(A)*B a1=randn(1); a2=rand(1); x=a1*Z(:,1)+a2*Z(:,2)+x0; norm(A*x-B) Z1=null(sym(A)) x01=sym(pinv(A)*B) syms a11 a22; xx=a11*Z(:,1)+a22*Z(:,2)+x01; 矩阵分解 1、LU分解和Cholesky分解 [L,U]=lu(A) D=chol(A) 2、奇异值分解 S=svd(A) 方程组的迭代解法 GS迭代 n = 10; A = rand(n)+ n*eye(n); b = A*ones(n,1); % For example, to solve Ax=b D = diag(diag(A)); L = -tril(A) + D; U = -triu(A) + D; DL= D-L; g = DL\b; tol = 1e-6; x = zeros(n,1); err = 1; while(errtol) xt = DL\(U*x) + g; err = norm(xt-x); %keyboard; x = xt; end x SOR迭代 n = 10; A = rand(n)+ n*eye(n); b = A*ones(n,1); % For example, to solve Ax=b D = diag(diag(A)); L = -tril(A) + D; U = -triu(A) + D; omega = 0.6; DL= D-omega*L; DU= (1-omega)*D+omega*U; g = omega* (DL\b); tol = 1e-6; x = zeros(n,1); err = 1; while(errtol) xt = DL\(DU*x) + g; err = norm(xt-x); %keyboard; x = xt; end * * 例 、 clear all 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)*B norm(A*x-B) x = -1.8000 2.4000 1.8667 -1.2667 3.8667 -3.2667 -2.1333 2.7333 ans = 7.4738e-015