线性方程组的直接法.ppt

1,第三章 解线性方程组的直接方法1 引言2 高斯消去法3 选主元素的高斯消去法4 矩阵的三角分解5 解三对角线方程组的追赶法6 解对称正定矩阵方程组的平方根法,2,1 引言,学习线性方程组数值解法的必要性科学计算中经常遇到线性方程组求解问题 如电路分析、分子结构、测量学、运筹学、流体力学、数值逼近及微分方程的数值解法等当线性方程组的阶数较大时，人工求解已不可能。当求解方法不得当时，即使计算机求解都很难实现。对20阶的线性方程组，用Cramer法则求解，乘除法的运算次数达到9.71020，若用每秒钟一亿次的计算机计算也要30万年；而用Gauss消去法求解，则只需乘除法次数为3060次，不需1秒钟就可计算出来。线性方程组求解方法的分类 直接法、间接法,3,矩阵的特征值和谱半径,4,特征值的性质及特征多项式,特征值的有关性质AT和A具有相同的特征向量及特征值；若A非奇异，则A1与A的特征值互为倒数，特征向量相同；相似矩阵具有相同的特征值。特征多项式,例题：P.22 例3.1 求矩阵的特征值及谱半径,5,对称正定矩阵,定义性质非奇异，且其逆矩阵也对称正定；所有特征值大于零；所有对角元也大于零；所有顺序主子式都大于零。,6,初等矩阵,定义性质,7,初等置换矩阵,初等置换矩阵的性质（1）是对称矩阵；（2）是正交矩阵；（3）行列式为1；（4）左乘矩阵A相当于将A的第i，j行互换，右乘A相当于将A的第i,j列互换。,8,初等下三角阵（Gauss变换矩阵）,定义：,性质：,9,Household变换（初等反射矩阵）,定义性质,10,2 高斯消去法 高斯消去法是一个古老的求解线性方程组的方法，但由它改进得到的选主元的高斯消去法则是目前计算机上常用的解低阶稠密矩阵方程组的有效方法。例1 用高斯消去法解方程组解 第1步：将方程（2.1）乘上（-3/2）加到方程（2.2）将方程（2.1）乘上（-1/2）加到方程（2.3）则得到与原方程等价的方程组,（2.1）,（2.3）,（2.2）,11,其中方程（2.4），（2.5）已消去了未知数。第2 步：方程（2.4）乘上2加到方程（2.5），消去（2.5）式中未知数，得到等价的三角形方程组由上述方程组，用回代的方法，即可求得原方程组的解。,（2.4）,（2.5）,12,若用矩阵来描述消去法的约化过程，即为这种求解过程，称为具有回代的高斯消去法。从上例看出，用高斯法解方程组的基本思想是用矩阵的初等变换将系数矩阵约化为具有简单形式的矩阵（上三角矩阵，单位矩阵等），从而容易求解。下面讨论求解一般线形方程组的高斯消去法，设有n个未知数的线性方程组：,13,引进记号,(2.7),(2.7)可用矩阵形式表示,(2.8),14,为了讨论方便，记假设 为非奇异矩阵（即设)。第1步（k=1）：设 计算乘数用 乘上（2.7）第一个方程，加到第i个中方程上去，即施行行初等变换：,15,(2.9),消去第2到第n个方程的未知数，得到（2.7）的等价方程组,其中（2.9）式中方框内元素为这一步需要计算的元素，计算公式为：,记为,16,第k步：继续上述消去过程，设第1步至第k-1步计算已 经完成，得到与原方程组等价的方程组。,（2.10）,17,记为,现进行第k步消元计算，设，计算乘数,用 乘（2.10）的第 k个方程加到第i个方程消去（2.10）中第i 个方程 的未知数 得到原方程组的等价方程组。,18,(2.11),简记为,其中 元素计算公式为：,19,重复上述约化过程，即 且设，共完成n-1步消元计算，得到与原方程组（2.7）等价的三角形方程组,(2.12),(2.13),20,用回代法，即可求得（2.13）的解，计算公式为：,（2.14）,元素 称为约化的主元素。将（2.7）约化为（2.13）的过程称为消元过程；（2.13）求解过程（2.14）称为回代过程，由消元过程和回代过程求解线性方程组的方法称为高斯消去法。,21,定理1（高斯消去法）设 其中 如果约 化的主元素 则可通过高 斯消去法（不进行交换两行的初等变换）将方 程组 约化为三角形矩阵方程组（2.13），且消元和求解公式为：1.消元计算,定理证明详见课本P.26定理2.2,22,2.回代计算,当A为非奇异矩阵时，也可能有某 但在第k列存在元素,23,于是可能通过交换（A，b）的第k行和第 行将 调到 位置，然后再进行消元计算。于是，在A为非奇异矩阵时，只要引进行交换，则高斯消去法可将原线性方程组 约化为三角形方程组（2.13），且通过回代法即可求得方程组的解。高斯消去法计算量:（1）消元计算：第k步 1.计算乘数：需要作（n-k）次除法运算；2.消元：需作 次乘法运算；3.计算：需作（n-k）次乘法运算；于是，完成全部消元计算共需作乘除运算 的次数为s：,24,（2）回代计算：共需要作n(n+1)/2次乘除运算。于是，用高斯消去法解 的计算量为共需作,（2.15）,次乘除运算。,25,下面比较用高斯消去法和用克莱姆（Cramer）法则解20阶方程组的计算量。,如果计算在每秒作10亿次乘除法运算的计算机上进行，那么用高斯消去法解20阶方程组约需要0.000003秒时间即可完成，而用克莱姆法则大约需 小时完成（大约相当于 年）。由此可知克莱姆法则完全不适用在计算机上求解高维方程组。,26,在计算机上用高斯消去法解低阶稠密矩阵线性方程组时要注意几点:(1)要用一个二维数组A(n,n)存放系数矩阵A的元素，用一维数组b(n)存放常数项b向量；(2)需要输入的数据：A，b，n(3)约化的中间结果 元素冲掉A元素，冲掉b，乘数 冲掉。,(4)在高斯消去法中一般要引进行交换；如果不存在 使,要输出方程没有唯一解 的信息。,27,例题：用Gauss消去法解方程组并求其系数矩阵行列式的值。,28,消去法与三角分解,Guass消元的过程从矩阵变换角度出发，实质上是进行了(n-1)次Gauss变换，即,29,3 选主元素的高斯消去法 用高斯消去法解 时，其中设A为非奇异矩阵，可能出现 情况，这时必须进行带行交换的高斯消去法。但在实际计算中即使 但其绝对值很小时，用 作除数，会导致中间结果矩阵 元素数量级严重增长和舍入误差的扩散,使得最后的计算结果不可靠。例2 设有方程组,30,解 精确解为方法1：用高斯消去法求解(用具有舍入的6位浮点数 进行运算),回代得到计算解。与精确解比较,这是一个很坏的结果。,31,回代求解。对于用具有舍入的6位浮点数进行运算，这是一个很好的计算结果。方法1计算失败的原因，是用了一个绝对值很小的数作除数，乘数很大，引起约化中间结果数量很严重增长，再舍入就使得结果不可靠了。,方法2 用具有行交换的高斯消去法(避免小主元)。,32,在采用高斯消去法解方程组时，小主元可能导致计算失败，故在消去法中应避免采用绝对值很小的主元素。对一般方程组，需要引进选主元的技巧，即在高斯消去法的每一步应该选取系数矩阵或消元后的低阶矩阵中选取绝对值最大的元素作为主元素，保持乘数 以便减少计算过程中舍入误差对计算解的影响；对同一个数值问题，用不同的计算方法，得到的结果的精度大不一样，一个计算方法，如果用此方法的计算过程中舍入误差得到控制，对计算结果影响较小，称此方法为数值稳定的；如果用此计算方法的计算过程中舍入误差增长迅速，计算结果受舍入误差影响较大，称此方法为数值不稳定。解数值问题时，应选择和使用数值稳定的计算方法，否则，如果使用数值不稳定的计算方法去解数值计算问题，就可能导致计算失败。,33,3.1完全主元素消去法 设有线性方程组 其中A为非奇异矩阵。方程组的增广矩阵为,第1步（k=1）:首先在A中选主元素,即选择，使,34,再交换A，b 的第1行与第 行元素，交换A的第1列与第 列元素(相当于交换未知数 与)，将 调到(1，1)位置(交换后为简单起见增广阵仍记为A，b，其元素仍记为)然后，进行消元计算。,第k步：继续上述过程，设已完成第1步到第k-1步 计算，A,b 约化为下述形式（仍记 元 素为 为元素为）：,35,第k步选主元区域,对于 做到（3）(1)选主元素：选取 使,36,（2）如果；则交换A,b第k行与第 行元素，若，则交换A的第k列与第 列元素。（3）消元计算（4）回代求解。经过上面的过程，即从第1步到n-1步完成选主元，交换两行，交换两列，消元计算，原方程组 约化为,37,其中 为未知数 调换后的次序。回代求解,3.2 列主元素消去法完全主元消去法是解低阶稠密矩阵方程组的有效方法，但完全主元素方法在选主元时要花费一定的计算时间；引入一种在实际计算中常用的部分选主元（列主元）消去法。列主元消去法即是每次选主元时，仅依次按列选取绝对值最大的元素作为主元素，且仅交换两行，再进行消元计算。,38,第k步选主元区域,设列主元素消去法已经完成第1步到第k-1步的按列选主元，交换两行，消元计算得到与原方程组等价的方程组,39,第k步计算如下：对于 做到（4）（1）按列选主元：即确定 使（2），则A为奇异阵，停止计算。（3），则交换A,b第 与第k 行元素。（4）消元计算：,40,（5）回代计算：,计算解 在常数项 内得到。例2的方法2就是列主元消去法。例3 用列主元素消去法解方程组,41,解 精确解为（舍入值）：,42,回代即得到计算解 本例是用具有舍入的4位浮点数进行运算，所得到的计算解还是比较准确的。完全主元素消去法框图（图61）,43,图6-1,完全主元高斯消去法方法框图,44,图6-2,列主元高斯消去法方法框图,45,图6-3,不选主元高斯消去法方法框图,46,用完全主元素法解，可用一整型数组 开始记录未知数 次序即，最后记录调整后未知数的足标。系数矩阵A存放在二维数组 内，常数项b存放在 内，解存放在数组 内。,47,4 矩阵的三角分解 现用矩阵理论来研究高斯消去法，设约化主元素 由于对A施行行的初等变换相当于用初等矩阵左乘于A，高斯消去法第1步：对应有其中,一、利用Gauss消元法实现矩阵的三角分解,48,第k步消元过程：对应有其中,k列,(4.1),49,利用递推公式（4.1）,则有,k列,(4.2),50,由（4.2）式得到,(4.3),其中,51,L为由乘数构成的单位下三角阵，U为上三角阵，（4.3）式表明,用矩阵理论来分析高斯消去法,得到一个重要结果,即在条件 下,高斯消去法实质上是将A分解为两个三角矩阵的乘积A=LU。二、矩阵三角分解实现的条件 显然，由本章高斯消去法及行列式性质知，如果,则有,52,反之，可用归纳法证明：如果A的顺序主子式则 总结上述讨论，得到下述重要定理。,其中,53,定理2（矩阵的三角分解）设，如果A的顺 序主子式 则A可分解为一个单位下三角阵与一个上三角 阵的乘积，即 且分解是唯一的。证明 仅就 来证明唯一性。设有 其中 为下三角角阵，为上三角阵。,（4.4）,54,上式右边为上三角阵，左边为单位下三角阵，故应为单位阵。即,由假设知 存在，于是从（4.4）可得,称矩阵的三角分解 为杜利特尔（Doolittle)分解。其中,55,在定理2条件下，同样可有三角分解其中L为下三角阵,U为单位上三角阵。称矩阵的这种分解为克劳特（Crout)分解。例4 由例1可得到A的三角分解,设有方程组。如果实现了,则求解问题化为：,56,求解两个三角形方程组是容易的。如果A为非奇异矩阵，则由列主元消去法知A通过行交换后可实现三角分解，即存有置换阵P使其中L为单位下三角阵，U为上三角阵。,即,57,三.矩阵的直接三角分解法 当矩阵A非奇异且顺序主子式，则A可做LU分解，即A=LU.这时方程组Ax=b可以转化为求解LUx=b，等价于以下两个三角形方程组 Lyb，及 Uxy.现在通过矩阵乘法直接求得L及U的元素。对比两边元素，直接可得,58,若己经求得U的第r-1行及L的第r-1列，则由矩阵乘法可得即从而可以求得U的第r行元素。L的第r列元素可由若,可得,59,求解三角形方程组的算法：对下三角方程组Ly=b,有 对上三角方程组Ux=y,有,60,例5 利用直接三角分解法求下列方程组的解,解法：将A进行三角分解,解方程组Ly=b，得y1=1，y2=-1，y3=4解方程组Ux=y，得x1=-1/2,x2=1,x3=0.,61,5 解三对角方程组的追赶法 在一些实际问题中，如用三次样条函数的插值问题，用差分解二阶线形常微分方程边值问题等，最后都导致解三对角方程组 即,（5.1),62,其中A满足条件,（5.2),对于具有条件（5.2)的方程组(5.1),可用下述的追赶法求解。追赶法具有计算量少，方法简单，算法稳定等特点。定理3 设三对角方程组，且A满足条件(5.2)，则A为非奇异矩阵。证明 用归纳法证明，显然，对n=2时有,63,现设定理对n-1阶满足条件（5.2)的三对角阵成立，求证对满足条件(5.2)的n阶三对角阵定理亦成立。因，消去法执行一步后得,显然,64,其中,于是由归纳法假设，则有,65,定理4 设，其中A为满足条件（5.2）的三对角阵，则A的所有顺序主子式都不为零，即,故,证明 由于A是满足（5.2）的n阵三对角阵，因此，A的任一个顺序主子阵亦是满足（5.2）的三对角阵，由定理3，则有：,66,即,于是，由矩阵的三角分解定理，则有,67,由矩阵乘法，可得计算待定系数 的计算公式，即,68,于是，得到解（5.1)的追赶法公式：(1)分解计算公式 A=LU,(2)求解 的递推公式,69,(3)求解 的递推公式,将计算 的过程称为追的过程,计算方程组解的过程称为赶的过程.追赶法解 仅需要5n-4次乘除运算。只需要用3个一维数组分别存贮A的系数 且还需要用三个一维数组保存计算的中间结果,70,例5 用追赶法解方程组,解(1)计算(2)计算(3)求解计算,71,6 解对称正定方程组的平方根法一、引言 在工程技术问题中,例如用有限元方法解结构力学中问题时,常常需要求解具有对称正定矩阵方程组,对于这种具有特殊性质系数矩阵,利用矩阵的三角分解法求解就得到解对称正定矩阵方程组的平方根法，这是一种解对称正定矩阵方程组的有效方法。二、对称正定矩阵 设有方程组,(6.1),其中,若 A满足下述条件：,72,（1）A对称：即（2）对任意非零向量 有 则称A为对称正定矩阵。对称正定矩阵A具有性质：a.设A为对称正定阵，则A的顺序主子式都大 于零，即,b.A的特征值满足：,73,三、平方根法 由设A为对称正定矩阵，则A有三角分解,(6.2),其中,74,75,由设,于是 为单位下三角阵，为上角阵，由矩阵三角分解的唯一性，则，从而对称正定矩阵A有唯一分解式由(6.2)式可知,(6.3),又因为，故,76,代入（6.3），则有,于是，对角阵D还可分解,其中 为下三角阵。,77,定理5（对称正定阵的三角分解）设A为n阶对称正定矩阵，则有三角分解（且 唯一）：（1），其中L为单位下三角阵，D 为对角阵，或（2），其中L为下三角阵且当限定L 的对角元素为正时这种分解是唯一的。这种矩阵分解称为乔来斯金（Cholesky）分解。分解计算 的递推公式及求解公式：,78,设有，其中 为对称正定阵，于是有三角分解,其中,79,由此得解对称正定矩阵方程组的平方根法计算公式：,由矩阵乘法，则有L的第一列元素为,同理，可确定L的第j列元素,80,（一）分解计算,（1）,（2）对于,（二）求解计算（3）求解,81,（4）求解,平方根法仅需计算L，计算量约为 次乘除法运算，大约为一般高斯消去法计算量的一半。,82,由分解公式有,（6.4）,（6.4）说明解 的平方根法所得的中间数据 是有界的，即 数量级不会增长。因此，虽然解对称正定矩阵方程组的平方根法没有进行选主元素，但平方根法是数值稳定的。因为A为对称矩阵，因此，电算时只需用一维数组存贮A的对角线以下元素，即,83,存区,用一维数组按行存贮：,84,解 精确解,（1）分解计算A=LLT,且矩阵A的元素 在一维数组中表示为：,计算得到的L的元素存放在A的对应位置。例6 用平方根法解方程组,85,（2）求解两个三角形方程组，求解 得到 求解，即得 的解,86,四、改进的平方根法,对称正定矩阵A可分解为A=LDLT其中,87,改进平方根法求解对称正定方程组,