《基于matlab的数值分析ppt培训课件》由会员分享，可在线阅读，更多相关《基于matlab的数值分析ppt培训课件（138页珍藏版）》请在金锄头文库上搜索。
1、基于MATLAB的数值分析,以软件MATLAB作为辅助工具介绍数值分析（科学与工程计算）的基本内容，注重讲授一些求解方程以及结果可视化的知识和技巧，使同学们能够有效地解决问题并处理计算结果。,内容包括： 1、MATLAB编程和绘图 2、数值分析的数学基础 3、数值算法在工程、科学中的应用,第一章 MATLAB入门,1、MATLAB的命令窗口工作都在此处完成。 2、怎样进行计算运算对象：矩阵算术运算符：+ - * .* ./ .倒除： ab=b/a变量与变量名：变量名和变量名类型不需声明。,3、数据显示格式默认格式：5位（format short)format long 16位format e
2、短的浮点格式format long e 长的浮点格式 4、清除命令clear:清除所有使用过的变量或某个（些）变量clc: 清除命令窗口 5、程序结构分支：if _else_end; if_elseif_end; if_break_end循环：for_end; while_end,6、读写输入数据：z=input(type youe input:) 键盘输入格式化输出：fprintf(e_format %12.5e n,vol) 7、数学函数 8、功能函数sort(x) sum(x) max(x) min(x) mod(x,y) rand(n) eval(s) 9、编程（编写M文件） 10、绘
3、图,第二章 数值代数,内容：数值代数就是研究有关矩阵计算的问题。主要包括：1、 线性代数方程组的求解；2、 矩阵特征值问题 要求：1、掌握用MATLAB求解的方法2、知道那些问题是困难的，那些问题是不可解的。,A=zeros(m,n) m行n列的零矩阵I=eye(n) n阶单位矩阵A=ones(m,n) 元素均为1A A的转置A(:,k) A(k,:) A(m1:m2,n1:n2)inv(A) A的逆size(A) A的大小hilb(n) Hilbert矩阵,2.1 矩阵,情况1：m=n(正规方程），最常见； 情况2：mn(超定方程）； 本节只介绍情况1。,MATLAB命令：,2.2 解线性代
4、数方程组的MATLAB命令,线性代数方程组并不总是数值可解的。只有当矩阵A的行列式不为零时才行！矩阵A的行列式即使不为零，但当很小或很大时，解的误差可能很大。,计算矩阵行列式的MATLAB命令：,2.4 病态问题有许多线性代数方程组理论上是可解的，但实际计算中由于受到舍入误差的影响而无法得到精确解。此类问题成为病态问题。病态问题的计算过程中，小的舍入误差或系数矩阵的微小变化都可能使解产生很大误差。（例子 P97),2.3 不可解问题,病态矩阵的一个重要标志是条件数：,MATLAB命令：,当矩阵是病态时，其条件数一定很大，但它并不能直接说明解的误差。线性方程组解的误差程度也取决于计算环境的精度。
5、条件数和行列式与计算环境是相互独立的。所以大条件数或小行列式未必意味无法直接精确求得线性方程组的解，它只意味着有很大误差可能。而实际上如果采用更高精度的计算环境则很可能得到非常满意的解。Hilbert矩阵是非常著名的病态矩阵（hilb(n),它经常用来检验算法的数值稳定性的好坏。,两种原因使我们想了解求解线性代数方程组的算法。一是实际工作中要用其它计算机语言（Fortran&C等）编写应用程序；二是MATLAB处理大型稀疏矩阵方程组显得很笨拙或无能为力。,2.5 线性代数方程组的求解方法（算法）,由线性代数的理论：,下面讨论实现过程：第一步：消元。进行到第k步时必有,a(k,k)作为主元，第k
6、行依次乘a(i,k)/a(k,k)加到第i行,i=k+1:n。 总共n-1步完成，k=1:n-1.,一、Gauss消去法,当a(k,k)=0，则上述消去法无法进行；或当其绝对值相对太小可能会出现大的计算误差。选主元法可避免这种情况。下面介绍常用的按列选主元的Gauss法。,列主元Gauss法,运算量 (只考虑乘除运算）： 第k步=n-k次除法+ n-k次乘法+（n-k)*(n-k)次乘法;,总的乘除运算量=,第二步：回代求解,%,二、LU分解法,LU分解的目的是将矩阵A转换为两个矩阵的乘积，即,好处是：对于线性方程组，如果需要多次求解不同的非齐次项，此时LU分解的效率将大大超过高斯消去法。,L
7、U分解的MATLAB命令：l,u=lu(A) 和l,u,p=lu(A),前面讲到的不选主元的高斯消去法和列主元高斯消去法将能实现LU分解。不选主元的高斯消去法用于下面两类矩阵肯定能成，即严格对角占优矩阵或对称正定矩阵，其他矩阵就难说了！列主元高斯消去法是解决一般中小型稠密矩阵方程组最有效的方法之一。下面讲解列主元高斯消去法实现LU分解的算法。,1、LU分解的代数理论,现在我们只要将列主元高斯消去法稍加改造即是LU分解的算法。列主元高斯消去法的矩阵表示：,2、LU分解算法,运算量和高斯消去法一样,3、三对角矩阵方程的追赶法,三对角且主对角严格占优矩阵方程是一类来源丰富的问题。比如，微分方程数值解
8、或样条插值等问题中的正规方程组。解这种问题必须考虑其矩阵稀疏的特征，减少算法的计算量。,三对角矩阵形如：,T的LU分解具有形式：,由T=LU推得：,最终解为f. 乘除运算量： 5n=O(n).,4、对称正定矩阵方程的cholesky分解法,对称正定矩阵方程的来源比较丰富，比如线性回归、拟合等问题。解这类问题必须考虑其矩阵对称的特征，减少算法的计算量。,由于A为对称正定矩阵，A必有cholesky分解：,三、线性代数方程组的迭代法,线性代数方程组的迭代法并不适用于所有问题，但它对一些特定类型的问题非常有效。当问题是大型稀疏矩阵方程时，高斯消去法的效率会变得非常低，而且有时还会超出内存要求。对于这
9、样的问题就需要使用迭代法。大系统问题的求解最终归结为大型稀疏矩阵方程。比如， 电网络、场方程的数值计算、运筹问题等。尽管迭代法的种类很多，这里只介绍其中的三种：1、Jacobi迭代；2、Gauss-Seidel迭代;3、超松弛迭代法（SOR)。,1、线性代数方程组的迭代法的一般理论,迭代算法简单，但问题是：1、能否保证算法收敛？2、能否充分利用矩阵的稀疏性，使运算量和存贮量尽量少。,定理一 迭代法收敛的充分必要条件是,定理二 迭代法收敛的充分必要条件是,定理三 迭代法收敛的充分条件是,收敛性定理,对于任意的初值,迭代矩阵B的构造,充分利用矩阵的稀疏性，使运算量和存贮量尽量少的办法就是要求迭代矩
10、阵B与原矩阵A有相同的稀疏结构。具体就是：,常用迭代法及其收敛性,定理4：SOR收敛的必要条件是02. 定理5：如果A是严格对角占优矩阵或不可分弱对角占优矩阵，则 1、Jacobi收敛； 2、gauss-seidel收敛； 3、当0=1 时，SOR必收敛。 定理6：如果A是对称正定矩阵，则 1、当2D-A正定时，Jacobi收敛； 2、gauss-seidel收敛； 3、当0eytat=r*r/(r*a*r);x=x+t*r;r=b-a*x; end,这一算法有简单易行，可充分利用A的稀疏性等特点，但当A的最大特征值远远大于最小特征值时收敛速度变得非常慢，以至于完全不适用。最速下降法并非最速！
11、下面的共轭梯度法可有效解决这一问题。,3、共轭梯度法,能否选出比负梯度方向更好的下降方向吗？能！,这一算法称为共轭梯度法，简称CG(conjugate gradient),实际计算中由于误差的影响，,之间的正交性很快就会,消失，以至于有限步内不能得到精确解，所以CG实际是一种 迭代算法。,CG是保稀疏且便于并行处理的算法，它是求解对称正定矩阵方程最优秀的算法之一。而且目前它已推广到一般矩阵方程，称为预优共轭梯度算法。,实用共轭梯度算法如下：,function x=cg(a,b,x,eyta) n=length(b); r=b-a*x; p=r; q0=r*r; while q0eytaw=a*
12、p;t=q0/(p*w);x=x+t*p;r=r-t*w;q=r*r;s=q/q0;p=r+s*p;q0=q; end,运算量：每一步迭代的乘除运算量为,2.6 矩阵特征值问题,一、幂法,MATLAB命令：D=eig(A)和X,D=eig(A),实用算法：,function lmta,x=mifa(a,epsl) m,n=size(a); x0=ones(n,1); while 1x=a*x0;lmta,m=max(abs(x);lmta=sign(x(m)*lmta;x=x/lmta;if abs(x-x0)epsl,break,end;x0=x; end,一、幂法,二、反幂法及其原点位移,
13、反幂法用来求A的按模最小的特征值。,实用算法：,function lmta,x=fanmifa(a,epsl) m,n=size(a); x0=ones(n,1); while 1x=ax0;lmta,m=max(abs(x);lmta=sign(x(m)*lmta;x=x/lmta;if max(abs(x-x0)epsl,break,end;x0=x; end lmta=1/lmta;,1、反幂法,2、带原点位移的反幂法,反幂法与“原点位移”相配合，求指定点附近的某个特征值和特征向量，并可用于加速幂法的收敛性。,另一方面，s离A的某个特征值越近，收敛越快。因此 不论用幂法求A的按模最大特征
14、值，还是利用反幂法求A的按模最小特征值，为了加快收敛，均可以用迭代m步后的近似值lmta作为最初始的位移值,实行动态位移迭代。,动态位移幂法,function lmta,x=dongtamifa(a,epsl0,epsl) lmta,x=mifa(a,epsl0);x0=x;lmta0=lmta; n=length(x0); while 1x=(a-lmta0*eye(n)x0;lmta,m=max(abs(x);lmta=sign(x(m)*lmta;x=x/lmta;lmta=1/lmta+lmta0;if max(abs(x-x0)epsl,break,end;x0=x;lmta0=lm
15、ta end,动态位移反幂法,function lmta,x=dongtaifanmifa(a,epsl0,epsl) lmta,x=fanmifa(a,epsl0);x0=x;lmta0=lmta; n=length(x0); while 1x=(a-lmta0*eye(n)x0;lmta,m=max(abs(x);lmta=sign(x(m)*lmta;x=x/lmta;lmta=1/lmta+lmta0;if max(abs(x-x0)epsl,break,end;x0=x;lmta0=lmta; end,三、幂法综述,1、幂法和反幂法只能用于求解可对角化矩阵的实数特征值和特征向量，不能求解复数特征值； 2、幂法和反幂法均为线性收敛，收敛速度由收敛因子决定，效率不高。 3、动态位移可以大大减小收敛因子加速收敛； 4、不适于求解全部特征值。 5、对称矩阵自然适用于幂法，此时采用2-范数标准化向量的算法至少平方收敛。,
