数值分析课件典型例题与习题2.ppt

三、四章内容提要三、四章内容提要典型例题分析典型例题分析数值分析典型例题典型例题 II化难为易化难为易化繁为简化繁为简化繁为简化繁为简初等行变换不改变方程组的解初等行变换不改变方程组的解1.交换矩阵的第交换矩阵的第i行与第行与第j行的位置行的位置2.以非零数以非零数k乘以矩阵的第乘以矩阵的第i行的行的每个元素每个元素3.把矩阵的第把矩阵的第i行的每个元素的行的每个元素的k倍倍加到第加到第j行的对应元素上去行的对应元素上去A(n 1)=Fn-1Fn-2F1 A其中其中Fk为为 Frobenius矩阵。矩阵。A=F1-1F2-1 Fn-1-1 A(n 1)直接方法直接方法:高斯消元法高斯消元法LU高斯消元法本质是矩阵的分解。矩阵分解矩阵分解(Top 10 Algorithms)(1)特征值分解特征值分解:A=CDC,C,D=eig(A)(3)LU分解分解:PA=LU,L,U,P=lu(A)(4)Cholesky分解分解:A=LLT,R=chol(A)(2)奇异值分解奇异值分解:A=USV,U,S,V=svd(A)(5)非负矩阵分解非负矩阵分解Learning the Parts of Objects by Non-negative Matrix Factorization,NatureNature,1999Demo1I=imread(monalisa.pgm);U,S,V=svd(double(I);s=diag(S);n1=5;Snew=diag(s(1:n1);zeros(size(s,1)-n1,1);figure,imshow(U*Snew*V,)n2=20;Snew=diag(s(1:n2);zeros(size(s,1)-n2,1);figure,imshow(U*Snew*V,)迭代法思想迭代法思想:Iterate:Tosayordoagainoragainandagain 迭代背后的思想是一种与传统思维模式截然不同的方式迭代背后的思想是一种与传统思维模式截然不同的方式,传统思维方式往往希望一遍做好传统思维方式往往希望一遍做好,一次成功一次成功;但是迭代开发意但是迭代开发意味着反复地做味着反复地做,不断地根据反馈进行调整。不断地根据反馈进行调整。迭代格式构造迭代格式构造收敛条件收敛条件(局部局部vs全局全局)中止准则中止准则加速加速(松弛思想松弛思想)Aitken加速方法加速方法 超松弛加速方法超松弛加速方法 共轭梯度法共轭梯度法的关键是构造一组两两共轭的关键是构造一组两两共轭的方向的方向(第第k步迭代生成共轭方向张成步迭代生成共轭方向张成k维子维子空间空间)。巧妙的是。巧妙的是共轭方向可以由上次搜索共轭方向可以由上次搜索方向和当前的梯度方向组合产生。方向和当前的梯度方向组合产生。现代迭代方法现代迭代方法(Top 10 Algorithms)最速下降法最速下降法思想简单思想简单,但收敛速度慢。本质但收敛速度慢。本质上因为负梯度方向函数下降快是局部性质。上因为负梯度方向函数下降快是局部性质。复杂性复杂性:高斯消元法共用乘法和除法次数为高斯消元法共用乘法和除法次数为n3/3+n2-n/3,常用记号常用记号O表示是多少阶的表示是多少阶的,则则O(n3/3)。注释注释2:复杂性对估计求解大型方程组所需的时间有用。例复杂性对估计求解大型方程组所需的时间有用。例如在一台特定的计算机上求解如在一台特定的计算机上求解n=500个方程的方程组所需个方程的方程组所需的时间我们可以通过求解一个的时间我们可以通过求解一个n=50个方程的方程组得到一个方程的方程组得到一个很好的猜测个很好的猜测,即对用掉的时间按比例放大即对用掉的时间按比例放大1000倍。倍。注释注释1:按按O的记法的记法,把把n的不同次幂相加的结果仅保留了最高的不同次幂相加的结果仅保留了最高次幂次幂,因为最高次幂决定了当因为最高次幂决定了当n趋近无穷时的极限形态。换趋近无穷时的极限形态。换而言之而言之,对于大的对于大的n,低阶项对算法的执行时间的估计没有太低阶项对算法的执行时间的估计没有太大影响大影响,仅需要近似估计执行时间时可以忽略不计。仅需要近似估计执行时间时可以忽略不计。常用的范数常用的范数:Hilbert矩阵条件数矩阵条件数:for i=1:10c(i)=cond(hilb(i),2);%vander(1:i)end,plot(1:10,c)范数范数(全局全局)问题的好与坏问题的好与坏 算法的快与慢算法的快与慢|X(k+1)X*|B|X(k)X*|范数的威力和魅力范数的威力和魅力:ekT=0 0 1 0 0=ImkekTFk1=I+mkekT(k=1,2,n1)F1-1F2-1 Fn-1-1=Fk1 Fj1=I+mkekT+mjejTkjF1-1F2-1 Fn-1-1=I+m1e1T+mn-1en-1T例例1.1.ekTmj=0,kj例例2.2.设设A为对称矩阵。为对称矩阵。高斯消元法一步后高斯消元法一步后,A约化为约化为 证明证明 B也是对称矩阵也是对称矩阵。约化主元不为零的判断约化主元不为零的判断定理定理3.1 约化主元约化主元 的充分必要条件是矩阵的充分必要条件是矩阵A各阶顺序主子式各阶顺序主子式不等于零。即不等于零。即例例4 4.严格对角占优矩阵的严格对角占优矩阵的约化主元约化主元ak,k(k-1)0(k=1,n)。例例3 3.严格主对角占优矩阵一定是非奇异的。严格主对角占优矩阵一定是非奇异的。严格对角占优矩阵严格对角占优矩阵:高斯消元法中约化主元高斯消元法中约化主元不等于零不等于零,Jacobi方法方法,GS方法和方法和SOR方法收敛。方法收敛。思考思考:若若A是严格对角占优矩阵是严格对角占优矩阵,当当0w=1时时,SOR方法收敛方法收敛。例例5 5.Ax=b b,其中其中A对称正定对称正定,问解此方程的问解此方程的Jacobi迭代法是否一定收敛？迭代法是否一定收敛？对称正定矩阵对称正定矩阵:直接法高斯消元法中约化直接法高斯消元法中约化主元大于零主元大于零,迭代法迭代法GS方法方法,SOR方法方法,最速下降最速下降方法和共轭梯度法收敛。方法和共轭梯度法收敛。例例6 6.例例7 7.例例8 8.病因病因是条件数大是条件数大,病症病症是什么呢是什么呢?例例9 9.矩阵的矩阵的Doolittle分解分解 A=LU,L为单位下三角矩阵为单位下三角矩阵,U为上三角矩阵。为上三角矩阵。a11=u11,a1n=u1na21=m21u11,an1=mn1u11对对A的元素的元素aij,当当jk和和 ik+1时时 m21u12+u22=a22,m21u1n+u2n=a2n u22=a22-m21u12,u2n=a2n-m21u1n m31u12+m32u22=a32,mn1u12+mn2u22=an2 m32=(a32-m31u12)/u22,mn2=(an2-mn1u12)/u22矩阵矩阵L和矩阵和矩阵U的元素计算公式的元素计算公式计算次序计算次序123456更新顺序更新顺序:先先行行后后列列列列除除行行不除不除旧元素减去旧元素减去所在所在行行和和列列前前k-1分量点乘的加和分量点乘的加和Doolittle分解分解：更新顺序更新顺序:先先列列后后行行行行除除列列不除不除旧元素减去旧元素减去所在所在行行和和列列前前k-1分量点乘的加和分量点乘的加和Crout分解分解：求矩阵的求矩阵的Doolittle分解分解三对角矩阵分解三对角矩阵分解步骤步骤I:三对角矩阵分解三对角矩阵分解 A=LU(k=2,3,n)下三角方程组下三角方程组 LY=f步骤步骤II:上三角方程组上三角方程组 UX=Y步骤步骤III:求解三对角方程组的追赶法求解三对角方程组的追赶法对称正定矩阵的对称正定矩阵的Cholesky分解分解(1)对于非零向量对于非零向量,xTAx总是正的总是正的;(2)A的顺序主子式全大于零的顺序主子式全大于零;(3)A的特征值全为正实数的特征值全为正实数。help chol对称正定矩阵对称正定矩阵:1)序列收敛序列收敛2)迭代矩阵谱半径小于迭代矩阵谱半径小于13)迭代矩阵特征值全小于迭代矩阵特征值全小于1 4)迭代矩阵范数小于迭代矩阵范数小于15)反证法反证法迭代法收敛证明思路迭代法收敛证明思路:例例1010.设设A是一个可逆矩阵是一个可逆矩阵,矩阵序列满足矩阵序列满足Xk+1=Xk(2IAXk)，（，（k=0,1,2,）则当则当 时时证明：由证明：由Xk+1=Xk(2IAXk)，得得IAXk+1=IAXk(2IAXk)=(IAXk)2 于是于是 IAXk=(IAXk-1)2 =(IAXk-2)22=例例1111.思路思路:(1)A1=B(I+R+R2+)；(2)任意给定任意给定n阶矩阵阶矩阵X0,由迭代格式由迭代格式Xk+1=XkR+B(k=0,1,2,)产生的矩阵序列产生的矩阵序列 Xk收敛到矩阵收敛到矩阵A-1；(3)对矩阵序列对矩阵序列 Xk,有误差估计式有误差估计式例例1212.设设A是是n阶可逆矩阵阶可逆矩阵,有有A的一个近似逆的一个近似逆B,令令 R=IAB。如果如果|R|q1,试证明试证明收敛的收敛的w取值范围。取值范围。例例1313.方程组方程组Ax=b,其中其中A是对称正定阵是对称正定阵,讨论迭代讨论迭代格式格式思路思路:例例14.14.思路思路:续例续例14.14.思路思路:例例15.15.思路思路:定理定理4.1 对对任意的任意的f和和任意的初始向量任意的初始向量X(0)迭代法迭代法 X(k+1)=BX(k)+f 收敛的收敛的充分必要充分必要条件是条件是例例16.16.思路思路:例例17.17.思路思路:-1/2a1/2收敛收敛,-1/2a1时系数矩阵正定。时系数矩阵正定。例例18.18.例例18.18.例例18.18.参考参考:Writting Fast Matlab Codes1.中小规模线性方程组中小规模线性方程组x=Ab;%Solves A*x=b2.(.(超超)大规模线性方程组大规模线性方程组(Preconditioned cg)作业作业题目题目1:从理论角度从理论角度(复杂度和收敛性复杂度和收敛性)比较各比较各 种方法的优劣。种方法的优劣。题目题目2:从数值角度比较各种方法的性能从数值角度比较各种方法的性能(公平公平)题目题目3:科研或生活中遇到的线性方程组？科研或生活中遇到的线性方程组？提示提示:参考代码参考代码:Matlab代码代码Iterative Solver1.如何生成方程组如何生成方程组A=gallery(poisson,n);b=ones(size(A,1),1);x0=zeros(size(A,1),1);2.如何生成方程组如何生成方程组矩阵市场矩阵市场:UF Sparse Matrix Collection:提示提示:3.更具挑战和趣味的例子更具挑战和趣味的例子图像编辑图像编辑 参考参考:Matlab代码代码Possion图像复原图像复原 参考参考:DeblurringImages:Matrices,Spectra,andFilteringGoogle PageRank 参考参考:NumericalComputingwithMATLAB数据挖掘和模式识别数据挖掘和模式识别 参考参考:MatrixMethodsinDataMiningandPatternRecognition高光谱图像解混高光谱图像解混 参考参考:SparseUnmixingofHyperspectralData医学成像医学成像 参考参考:MathematicalMethodsinImageReconstruction迭代法思想迭代法思想:Tostarttoday,tostartnowandtoship,makeitahabit 迭代背后的思想是一种与传统思维模式截迭代背后的思想是一种与传统思维模式截然不同的方式然不同的方式,传统思维方式往往希望一遍做好传统思维方式往往希望一遍做好,一次成功一次成功;但是迭代开发意味着反复地做但是迭代开发意味着反复地做,不断不断地根据反馈进行调整。地根据反馈进行调整。