数学建模线性方程组的数值解法.pptVIP

* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 线性方程组数值解法的MATLAB实现 若A为可逆方阵，输出原方程的解x 若A为n?m矩阵（nm），且ATA可逆，输出原方程的最小二乘解x 线性方程组数值解法的MATLAB实现 [x,y,p]=lu(A) 若A可逆，输出x为单位下三角阵L， y为上三角阵U，p为一交换阵P，使PA=LU. u =chol(A) 对正定对称矩阵A的Cholesky分解，输出u为上三角阵U，使A=UTU 2. 矩阵LU分解 [x,y]=lu(A) 若A可逆且顺序主子式不为零， 输出x为单位下三角阵L，y为上三角阵U，使A=LU； 若A可逆，x为一交换阵与单位下三角阵之积. 例. 解 A=[10 3 1;2 -10 3;1 3 10], b=[14 -5 14], x=A\b, [L1,U1]=lu(A); L1,U1, A1=L1*U1, [L2,U2,P]=lu(A); L2,U2,P, A2=L2*U2, A3=inv(P)*A2 并对系数矩阵作LU分解 shiyan51 若第1个方程改为 3x2+x3=14 结果如何 当n很大时Hilbert矩阵呈病态 线性方程组数值解法的MATLAB实现 3. 范数 条件数 特征值 n=norm(x) 输入x为向量或矩阵，输出为x的2-范数 c=cond(x) 输入x为矩阵, 输出为x的2-条件数 r=rcond(x) 输入x为方阵, 输出为x条件数倒数 e=eig(x) 输入x为矩阵，输出x的全部特征值 H=hilb(5), h=rats(H), b=ones(5,1); x=H\b; b(5)=1.1; x1=H\b; [x,x1], n1=cond(H), n2=rcond(H), 例：观察Hilbert矩阵的病态性 例. Hx=b, 其中 H=hilb(5), b=[1,…1]T shiyan52 x x1 1.0e+003 * 0.0050 0.0680 -0.1200 -1.3800 0.6300 6.3000 -1.1200 -9.9400 0.6300 5.0400 cond(H)=4.7661e+005 1. 提取（产生）对角阵 v=diag(x) 输入向量x，输出v是以x为对角元素的对角阵；输入矩阵x，输出v是x的对角元素构成的向量； 例：v=diag(diag(x)) 输入矩阵x，输出v是x的对角元素构成的对角阵，可用于迭代法中从A中提取D。 2. 提取上（下）三角阵 其他相关的MATLAB函数 y=triu(x) 输入矩阵x，输出v是x的上三角阵； v=tril(x) 输入矩阵x，输出v是x的下三角阵； v=triu(x,1) 同上，但对角元素为0，可从A中提取U; v=tril(x,-1) 同上，但对角元素为0，可从A中提取L。 例. 用迭代法解 shiyan53 MATLAB对稀疏矩阵的处理: 进行大规模计算的优点 a=sparse(r,c,v,m,n) 在第r行、第c列输入数值v，矩阵共m行n列，输出a为稀疏矩阵，只给出(r,c)及v aa=full(a) 输入稀疏矩阵a，输出aa为满矩阵（包含零元素） a=sparse(2,2:3,8,2,4), aa=full(a), a =(2,2) 8 aa= 0 0 0 0 (2,3) 8 0 8 8 0 输出 n=500;b=[1:n]; a1=sparse(1:n,1:n,4,n,n); a2=sparse(2:n,1:n-1,1,n,n); a=a1+a2+a2; tic;x=a\b;t1=toc aa=full(a); tic;xx=aa\b;t2=toc y=sum(x) yy=sum(xx) 例. 分别用稀疏矩阵和满矩阵求解Ax=b, 比较计算时间 设 0 0 t1, t2相差巨大，说明用稀疏矩阵计算的优点 （y=yy 用于简单地验证两种方法结果的一致） shiyan54 用MATLAB作线性最小二乘拟合 1. 作多项式 f(x)=a1xm+ …+amx+am+1拟合,可利用已有程序: a=polyfit(x,y,m) 输入:数据x,y (同长度数组);m (拟合多项式次数) 输出:系数a=[a1, …am , am+1] (数组)。 2. 对超定方程组 仍用 可得最小二乘意义下的解 多项式