科学计算与数学建模第五章.ppt

《科学计算与数学建模第五章.ppt》由会员分享，可在线阅读，更多相关《科学计算与数学建模第五章.ppt（82页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、科学计算与数学建模科学计算与数学建模科学计算与数学建模科学计算与数学建模中南大学数学科学与计算技术学院中南大学数学科学与计算技术学院 小行星轨道方程计算问题小行星轨道方程计算问题第五章第五章 小行星轨道方程计算问题小行星轨道方程计算问题 线性方程组求解的直接法线性方程组求解的直接法小行星轨道方程问题小行星轨道方程问题 5.1线性方程组直接解法概述线性方程组直接解法概述5.2直接解法直接解法5.3小行星轨道方程问题的模型求解小行星轨道方程问题的模型求解5.45.1.1 5.1.1 问题的引入问题的引入v一天文学家要确定一颗小行星绕太阳运行的轨道，他在轨道平面内一天文学家要确定一颗小行星绕太阳运行

2、的轨道，他在轨道平面内建立以太阳为原点的直角坐标系，其单位为天文测量单位，在建立以太阳为原点的直角坐标系，其单位为天文测量单位，在5 5个个不同的时间对小行星作了不同的时间对小行星作了5 5次观察，测得轨道上的次观察，测得轨道上的5 5个点的坐标数据个点的坐标数据如下表：如下表：表表 5.1.1 5.1.1 轨道上的轨道上的5 5个点的坐标数据个点的坐标数据试确立小行星的轨道方程，并画出小行星的运动轨线图形。试确立小行星的轨道方程，并画出小行星的运动轨线图形。123455.7646.2866.7597.1687.4080.6481.2021.8232.5263.360 小行星轨道方程问题小行星

3、轨道方程问题 5.1.2 5.1.2 模型的分析模型的分析 v由开普勒第一定律知，小行星轨道为一椭圆，设椭圆的一般方程为：由开普勒第一定律知，小行星轨道为一椭圆，设椭圆的一般方程为：，需要确定系数，需要确定系数 利用已知的数据，不妨设利用已知的数据，不妨设 欲确定系数欲确定系数 等价于求解一个线性方程组：等价于求解一个线性方程组：可写成矩阵的形式：可写成矩阵的形式：v其中，其中，5.1.3 5.1.3 模型的假设模型的假设 v假设：假设：(1)(1)小行星轨道方程满足开普勒第一定律；小行星轨道方程满足开普勒第一定律；(2)(2)以上所测得数据真实有效。以上所测得数据真实有效。5.1.3 5.1

4、.3 模型的建立模型的建立 v该问题的模型为：该问题的模型为：可见，解答上述问题就是对线性方程进行求解可见，解答上述问题就是对线性方程进行求解。5.2 5.2 线性方程组直接解法概述线性方程组直接解法概述v 直接法：直接法：利用一系列递推公式计算有限步，能直接得到方程组的精利用一系列递推公式计算有限步，能直接得到方程组的精确解。当然，实际计算结果仍有误差，譬如舍入误差。舍入误差的积累确解。当然，实际计算结果仍有误差，譬如舍入误差。舍入误差的积累有时甚至会严重影响解的精度有时甚至会严重影响解的精度v 求解线性方程组最基本的一种直接法是求解线性方程组最基本的一种直接法是消去法消去法。这是一个众所周

5、知。这是一个众所周知的古老方法，但用在现代电子计算机上仍然十分有效。消去法的基本思的古老方法，但用在现代电子计算机上仍然十分有效。消去法的基本思想是，通过将一个方程乘或除以某个常数，以及将两个方程思想是，通想是，通过将一个方程乘或除以某个常数，以及将两个方程思想是，通过将一个方程乘或除以某个常数，以及将两个方程相加减这两种手续，过将一个方程乘或除以某个常数，以及将两个方程相加减这两种手续，逐步减少方程中的变元的数目，最终使每个方程仅含一个变元，从而得逐步减少方程中的变元的数目，最终使每个方程仅含一个变元，从而得出所求的解。其中出所求的解。其中高斯消去法高斯消去法是广泛应用的方法，其求解过程分为

6、消元是广泛应用的方法，其求解过程分为消元过程和回代过程两个环节。消元过程将所给的方程组加工成上三角方程过程和回代过程两个环节。消元过程将所给的方程组加工成上三角方程组。所归结的方程组再通过回代过程得出它的解。高斯消去法由于添加组。所归结的方程组再通过回代过程得出它的解。高斯消去法由于添加了回代的过程，算法结构稍复杂，但这种算法的改进明显减少了计算量。了回代的过程，算法结构稍复杂，但这种算法的改进明显减少了计算量。v 直接法比较适用于中小型方程组。对高阶方程组，即使系数矩阵是直接法比较适用于中小型方程组。对高阶方程组，即使系数矩阵是稀疏的，但在运算中很难保持稀疏性，因而有存储量大，程序复杂等不稀

7、疏的，但在运算中很难保持稀疏性，因而有存储量大，程序复杂等不足。足。5.3.1 5.3.1 高斯消去法高斯消去法v 消去法是一个古老的求解线性方程组的方法。由它改进得消去法是一个古老的求解线性方程组的方法。由它改进得选主元法选主元法是目前计算机上常用的有效的求解低阶稠密矩阵线性方程组是目前计算机上常用的有效的求解低阶稠密矩阵线性方程组的方法。的方法。用用 消去法解方程组消去法解方程组 直接解法直接解法v解：第解：第 1 1 步，步，加到加到 ，加到加到 ，得等价方程组：得等价方程组：第第 2 2 步，步，加到加到 得等价的方程组：得等价的方程组：第第 3 3 步，回代法求解步，回代法求解 即可

8、求得该方程组的解为：即可求得该方程组的解为：用矩阵法描述的约化过程即为：用矩阵法描述的约化过程即为：v 这种求解过程称为具有回代的这种求解过程称为具有回代的高斯消去法高斯消去法。v 此此例例可可见见GaussGauss消消去去法法的的基基本本思思想想是是：用用矩矩阵阵A A的的初初等等行行变变换换将将系系数数矩矩阵阵化化为为具具有有简简单单形形式式的的矩矩阵阵（如如上上三三角角阵阵，单单位位矩矩阵阵等等），而而三三角角形形方方程组是很容易回代求解的。程组是很容易回代求解的。一般的，设有个未知数的线性方程组为：一般的，设有个未知数的线性方程组为：设设 则则 化为：化为：为方便，为方便，则消去法为

9、：则消去法为：第第1 1步：步：计算计算 用用 乘乘 的第一的第一方程加到第方程加到第 个方程中去个方程中去 （即实行行的初等变换）（即实行行的初等变换），消去第，消去第2 2个到第个到第n n个方程中的未知数个方程中的未知数 得与得与 等价方程组：等价方程组：记为：记为：其中其中 式中元素式中元素 为进一步需要计算的元为进一步需要计算的元素，公式为：素，公式为：v第第 ，步，继续上述过程消元。设第步，继续上述过程消元。设第 步到第步到第 步计步计算已完成，得到与原方程组等价的方程组：算已完成，得到与原方程组等价的方程组：记为记为 ，下面进行第，下面进行第 步消元法：步消元法：v设设 ，计算乘

10、数，计算乘数 v用用 中第中第 个方程加到第个方程加到第 个方程个方程 消去消去 中第中第 个方程个方程 的未知数的未知数 得到与原方程组等价的方程组：得到与原方程组等价的方程组：v记为记为 其中其中 中元素计算公式为：中元素计算公式为：最后，重复上述过程，即最后，重复上述过程，即 且设且设 共完成共完成 步消元计算，得到与步消元计算，得到与 等价的三角形方程组。等价的三角形方程组。再用回代法求解再用回代法求解 的解，计算公式为：的解，计算公式为：元素元素 称为称为约化的主元素约化的主元素。将。将 化为化为 的过程称的过程称为为消元过程消元过程。由消元过程和回代过程求解线性方程组的方法称为由消

11、元过程和回代过程求解线性方程组的方法称为Gauss Gauss 消去法消去法。的求解过程的求解过程 称为称为回代过程回代过程。v定理（定理（GaussGauss消去法）消去法）设设 若约化的主元素若约化的主元素 则可通过则可通过GaussGauss消去法（不进行两行的初等变换消去法（不进行两行的初等变换两行交换位置）将方程两行交换位置）将方程组化为等价的三角形方程组组化为等价的三角形方程组 。消元和求解的计算公式为：。消元和求解的计算公式为：1 1、消元计算、消元计算 2 2、回代计算、回代计算5.3.2 5.3.2 矩阵的三角分解矩阵的三角分解v 下面用矩阵理论进一步来分析下面用矩阵理论进一

12、步来分析 GaussGauss消去法，设约化主元素消去法，设约化主元素 由于对由于对 实行的初等变换相当于用初等矩阵左乘实行的初等变换相当于用初等矩阵左乘于是，于是，消去法第消去法第1 1步：步：,则有：则有：v v其中：其中：v （为初等三角矩阵）为初等三角矩阵）vGauss消去法第消去法第k步消元过程：步消元过程：v则有则有v其中：其中：利用递推公式则有：利用递推公式则有：由由 得得:v其中其中vL L为由乘数构成的下三角阵为由乘数构成的下三角阵,U,U为上三角矩阵，为上三角矩阵，表明，表明，用矩阵理论来分析用矩阵理论来分析GaussGauss消去法，得到一个重要结果，即消去法，得到一个重

13、要结果，即在在 条件下条件下GaussGauss消去法实质上是消去法实质上是A A将分将分解成两个三角矩阵的解成两个三角矩阵的v显然，可由显然，可由GaussGauss消去法及行列式性质可知，如果消去法及行列式性质可知，如果 v则有则有 其中其中 为顺序主为顺序主子式子式v反之，可用归纳法证明：如果反之，可用归纳法证明：如果A A的顺序主子式的顺序主子式 满足：满足：v则则 v定理定理 5.3.2 5.3.2（矩阵的三角分解）（矩阵的三角分解）设设 ，如果，如果 的顺序的顺序主子式有主子式有 ，则，则 可分解为一个单位下三角矩阵可分解为一个单位下三角矩阵与一个上三角矩阵的乘积，即与一个上三角矩

14、阵的乘积，即 且分解是唯一的且分解是唯一的。v证明证明 现仅就现仅就 来证明唯一性，存在性上面已证。来证明唯一性，存在性上面已证。假若假若 且对且对 非奇异时考虑，非奇异时考虑，为单位下三为单位下三角阵，角阵，为上三角阵，由假设知为上三角阵，由假设知 存在（因为存在（因为 可逆可逆 故故 可逆），从而由可逆），从而由 有有 ，上式右端，上式右端为上三角阵，左边为单位下三角阵，因此左右两端应为单位矩阵。故为上三角阵，左边为单位下三角阵，因此左右两端应为单位矩阵。故 即分解是唯一的。即分解是唯一的。v 称矩阵的三角分解称矩阵的三角分解 （杜利特尔）分解。（杜利特尔）分解。其中其中v在以上定理条件下

15、，同样可有下面的三角分解：在以上定理条件下，同样可有下面的三角分解：，其中，其中L L为为下三角矩阵，下三角矩阵，U U为单位上三角矩阵，称之为为单位上三角矩阵，称之为CroutCrout（克劳特）分解（克劳特）分解。v如前例中系数矩阵如前例中系数矩阵A A的分解为：的分解为：现设现设 ，若如分解，若如分解 ，则，则 而求解这两个三角形方程组是很容易的。而求解这两个三角形方程组是很容易的。v v 设设A A为为n n阶非奇异矩阵，则用阶非奇异矩阵，则用GaussGauss消去法解消去法解 所需要的乘除所需要的乘除法次数及加减法的次数分别为：法次数及加减法的次数分别为：Gauss Gauss消去

16、法的计算量消去法的计算量v但如果用但如果用 （克莱姆）法则解（克莱姆）法则解 ，就需要计算，就需要计算 个个n n阶阶行列式，若行列式用子式展开，总共需要行列式，若行列式用子式展开，总共需要 次乘法，如次乘法，如 时时 消去法需要消去法需要430430乘除法，而克莱姆法则却需要乘除法，而克莱姆法则却需要3991680039916800次乘次乘法，由此可见，法，由此可见，法则方程组的工作量太大，不便于使用。如果法则方程组的工作量太大，不便于使用。如果计算是在每秒作计算是在每秒作 次乘除法计算机上进的，那么用次乘除法计算机上进的，那么用 消去法解消去法解2020阶方程组约需秒即可完成，而用阶方程组

17、约需秒即可完成，而用 法则大约需法则大约需 小时才能完成（大约相当于小时才能完成（大约相当于 年）可见，年）可见，法则完全不适于在法则完全不适于在计算机上求解高维方程组。计算机上求解高维方程组。v用用GaussGauss消去法解消去法解 时，设时，设A A非奇异，但可能出现非奇异，但可能出现 ，这时必须进行行交换的，这时必须进行行交换的GaussGauss消去法，但在实际计算中消去法，但在实际计算中即使即使 ，但其绝对值很小时，用，但其绝对值很小时，用 作除数，会导致作除数，会导致中间结果矩阵中间结果矩阵 的元素数量级严重增长和舍入误差的扩的元素数量级严重增长和舍入误差的扩散，使最后结果不可靠

18、。散，使最后结果不可靠。设有方程组设有方程组5.3.4 Gauss主元素消去法主元素消去法v解解 方程组得解方程组得解v方法一：用方法一：用Gauss消去法求解。消去法求解。（用具有舍入的（用具有舍入的3位浮点数进行运算）位浮点数进行运算）回代得解回代得解 ，与精确解，与精确解比较，这结果很差。比较，这结果很差。v方法二：用具有行交换的方法二：用具有行交换的Gauss消去法（避免小消去法（避免小主元）主元）回代得解：回代得解：v这个解对于具有舍入的这个解对于具有舍入的3位浮点数进行计算，是一位浮点数进行计算，是一个很好的结果。个很好的结果。v方法一计算失败的原因，是用了一个绝对值很小方法一计算

19、失败的原因，是用了一个绝对值很小的数作除数，乘数很大，引起约化中间结果数量的数作除数，乘数很大，引起约化中间结果数量误差很严重增长，再舍入就使得计算结果不可靠误差很严重增长，再舍入就使得计算结果不可靠了。了。v这个例子告诉我们，在采用高斯消去法解方程组这个例子告诉我们，在采用高斯消去法解方程组时，小主元可能导致计算失败，故在消去法中应时，小主元可能导致计算失败，故在消去法中应避免采用绝对值很小的主元素。避免采用绝对值很小的主元素。v对一般矩阵方程组，需要引进主元的技巧，即在对一般矩阵方程组，需要引进主元的技巧，即在高斯消去法的每一步应该选取系数矩阵或消元后高斯消去法的每一步应该选取系数矩阵或消

20、元后的低阶矩阵中的绝对值最大的元素作为主元素，的低阶矩阵中的绝对值最大的元素作为主元素，保持乘数保持乘数 以便减少计算过程中以便减少计算过程中的舍入误差对计算解的影响。的舍入误差对计算解的影响。v 这个例子还告诉我们，对同一个数值问题，这个例子还告诉我们，对同一个数值问题，用不同的计算方法，得到的精度大不一样，一个用不同的计算方法，得到的精度大不一样，一个计算方法，如果用此方法的计算过程中舍入误差计算方法，如果用此方法的计算过程中舍入误差得到控制，对计算结果影响较小，称此方法为数得到控制，对计算结果影响较小，称此方法为数值稳定的；否则，如果用此计算方法的计算过程值稳定的；否则，如果用此计算方法

21、的计算过程中舍入误差增长迅速，计算结果受舍入误差影响中舍入误差增长迅速，计算结果受舍入误差影响较大，称此方法为数值不稳定。因此，我们解数较大，称此方法为数值不稳定。因此，我们解数值问题时，应选择和使用数值稳定的计算方法，值问题时，应选择和使用数值稳定的计算方法，否则，如果使用数值不稳定的计算方法去解数值否则，如果使用数值不稳定的计算方法去解数值计算问题，就可能导致计算失败。计算问题，就可能导致计算失败。v设有线性方程组，其中为非奇异矩阵。方程组得增广矩阵为设有线性方程组，其中为非奇异矩阵。方程组得增广矩阵为 第第 步步 ：首先在：首先在 中选主元素，即选择中选主元素，即选择 使使 ，再交换，再

22、交换 的第一行与第的第一行与第 行元素，交换行元素，交换 的第一列与第的第一列与第 列元素（相当于交换未知数列元素（相当于交换未知数 与与 ），将），将 调到调到 位置位置（交换后的增广矩阵为（交换后的增广矩阵为 ，其元素仍记为，其元素仍记为 ，然，然后进行消元法计算。后进行消元法计算。5.3.5 完全主元素消去法完全主元素消去法v第第 步：继续上述过程，设已完成第步：继续上述过程，设已完成第 步到第步到第 步计算，步计算，约化为下述形式（为简单起见，仍记约化为下述形式（为简单起见，仍记 元素为元素为 元素元素为为 ）：）：域 v于是第于是第 步计算：对于步计算：对于 按下述步骤从按下述步骤从

23、 计算到计算到第第 步主元区域步主元区域（1）选主元素：选择）选主元素：选择 使使 （2）如果）如果 ，则交换，则交换 第第 行与第行与第 行元素，如果行元素，如果 ，则交换，则交换 的第的第 列与第列与第 列元素。列元素。（3）消元计算）消元计算（4）回代求解。经过上面的过程，即从第）回代求解。经过上面的过程，即从第 步到第步到第 步完成选主元，交换两行，交换两列，消元计算，步完成选主元，交换两行，交换两列，消元计算，原方程组约化为：原方程组约化为：其中其中 为未知数为未知数 调换调换后的顺序。回代求解：后的顺序。回代求解：列主元消去法列主元消去法v 完全主元素消去法完全主元素消去法是解低阶

24、稠密矩阵方程组的有效是解低阶稠密矩阵方程组的有效方法，但完全主元素方法在选主元时要花费一定的时间。方法，但完全主元素方法在选主元时要花费一定的时间。现介绍一种在实际计算中常用的部分选主元，（即列主现介绍一种在实际计算中常用的部分选主元，（即列主元）消去法。元）消去法。列主元消去法列主元消去法，即是每次选主元时，仅依，即是每次选主元时，仅依次按列选取绝对值最大的元素作为主元素，且仅交换两次按列选取绝对值最大的元素作为主元素，且仅交换两行，再进行消元计算。行，再进行消元计算。设列主元消去法已经完成第设列主元消去法已经完成第 步到第步到第 步的按步的按列选主元，交换两行，消元计算得到与原方程组等价的

25、列选主元，交换两行，消元计算得到与原方程组等价的方程组：方程组：v 第第 步步选选主主元元区区域域v第第 步计算如下：对于步计算如下：对于 ,按按下述步骤从下述步骤从 计算到计算到 按列主元，即确定按列主元，即确定 使使 如果如果 ，则，则 为非奇异，停止计算。为非奇异，停止计算。如果如果 ，则交换，则交换 第第 行第行第 行行元素元素 消元计算：消元计算：（5）回代计算：回代计算：计算解计算解 在常数项内在常数项内 得到得到 用列主元消去法求解方程组用列主元消去法求解方程组解解 精确解为（舍入值）：精确解为（舍入值）：回代计算得到计算解：回代计算得到计算解：本例是具有舍入的位浮点数进行运算，

26、所得的计算解还是比较准本例是具有舍入的位浮点数进行运算，所得的计算解还是比较准确的。确的。下面是完全主元素消去法框图（图下面是完全主元素消去法框图（图 ）输入输入 选主元素选主元素 否否是是否否交换行交换行输出输出 否否交换列且交换列且消元计算消元计算回代回代（当 ）调整求解调整求解输入计算解输入计算解图图5.3.1 完全主元素消去法框图完全主元素消去法框图v用完全主元素消元法解用完全主元素消元法解 ，可用一整型数组，可用一整型数组 开始记录未知数开始记录未知数 次序，即次序，即 ，最后记录调整后未知数的足标。系数阵，最后记录调整后未知数的足标。系数阵 存在二维存在二维数组数组 内，常数项内，

27、常数项 存在存在 内，解内，解存在数组存在数组 内。内。v例例 5.3.4 若在计算过程中，只取若在计算过程中，只取3位有效数字，试用列位有效数字，试用列主元素法求解：主元素法求解：v解解 第一步，选第一步，选 为主列元，将为主列元，将 对对调位置，调位置，v第二步，选第二步，选 为列主元，不需换行，为列主元，不需换行，v由由 回代求解得：回代求解得：与原方程组准确解与原方程组准确解 比较，比较，可知，本题用可知，本题用3位有效数字计算的列元素法是相当精确的。大量实位有效数字计算的列元素法是相当精确的。大量实践表明：列主元法为解线性方程组的精确方法。践表明：列主元法为解线性方程组的精确方法。v

28、下面用矩阵运算来描述列主元素法下面用矩阵运算来描述列主元素法v记记 是初等排列阵（由交换单位矩阵是初等排列阵（由交换单位矩阵 的第的第 行与行与 行所得）。则列主元素法为：行所得）。则列主元素法为：v其中其中 的元素满足的元素满足 ，由，由 得：得：简记为：简记为：，其中，其中 。下面考察下面考察 时的时的 其中，其中，则由排列阵性质（左乘矩阵是对矩阵进行行变化。）则由排列阵性质（左乘矩阵是对矩阵进行行变化。）v 已知已知 为单位下三角阵，其元素绝对值为单位下三角阵，其元素绝对值 ，记，记 。由。由 得：得：。其中。其中 为排列阵，为排列阵，为为单位下三角阵，单位下三角阵，为上三角阵。为上三角

29、阵。v这表明，对这表明，对 应用列主元素法相当于对应用列主元素法相当于对 先进行先进行一系列行变换后对一系列行变换后对 再应用再应用 消去法。再实消去法。再实际计算中我们只能在计算过程中，做关于行的变换。有结论：际计算中我们只能在计算过程中，做关于行的变换。有结论：v定理定理5.3.4 （列主元素三角分解定理）（列主元素三角分解定理）v 若若 为非奇异性矩阵，则存在排列矩阵为非奇异性矩阵，则存在排列矩阵 使使 。其。其中中 为单位下三角阵，为单位下三角阵，为上三角阵。为上三角阵。v 存放在存放在 的下三角部分，的下三角部分，存放在存放在 的上三角部分。由的上三角部分。由整数型数组记录可知整数型

30、数组记录可知 的情况。的情况。5.3.7 Gauss-Jordan 消去法消去法v 消去法总是消去对角线下方的元素。现考虑一种消去法总是消去对角线下方的元素。现考虑一种修正，即消去对角线下方和上方的元素。这即为修正，即消去对角线下方和上方的元素。这即为 消去法。消去法。设用设用 消去法已完成消去法已完成 步，于是步，于是 化为化为等价方程组等价方程组 其中：其中：v在第在第 步计算时，考虑对上述矩阵的第步计算时，考虑对上述矩阵的第 行上、下都进行消元行上、下都进行消元计算计算1、按列选主元素，即定义按列选主元素，即定义 使使2、换行（当换行（当ikk）。交换）。交换A,b第第k行与第行与第ik

31、行元素。行元素。3、计算乘数计算乘数 (可存放在可存放在 的单元中的单元中)（）4、消元计算计算主元素消元计算计算主元素 v5、计算主元素计算主元素 上述过程全部执行完后有：上述过程全部执行完后有：这表明这表明 用方法将用方法将 约化为单位矩阵，计算解就在常数项约化为单位矩阵，计算解就在常数项位置得到，因此用不着回代求解。用位置得到，因此用不着回代求解。用 方法接方程组的计算方法接方程组的计算量大约需要量大约需要 次乘除法，要比次乘除法，要比 消去法大些。但用消去法大些。但用 方法求一个矩阵的逆矩阵还是比较合适的。方法求一个矩阵的逆矩阵还是比较合适的。v定理定理（G-J法求逆矩阵）法求逆矩阵）

32、设设 的逆矩阵的逆矩阵 ，即，即求求 阶矩阵阶矩阵 ，使，使 其中其中 为单位矩阵，将为单位矩阵，将 按列写成按列写成 ，为为列向量，列向量，为单位列向量。于是求解为单位列向量。于是求解 ，等，等价于求解价于求解 个方程组：个方程组：，所以我们可以用所以我们可以用 法求解法求解 。v 对对 求求 。v解解 故：故：直接三角分解直接三角分解v为求解为求解 将将 进行进行 分解。即分解。即 将原问题转化为求解两个三角形方程组将原问题转化为求解两个三角形方程组 求求 求求 。v 不选主元的三角分解法不选主元的三角分解法v设设 且有且有 其中其中 为单位下三角为单位下三角阵，阵，上三角阵。即上三角阵。

33、即 v下面说明：下面说明：的元素可由的元素可由 步直接计算出来，其中步直接计算出来，其中第第 步定出步定出 的第的第 行和行和 的第的第 列元素。列元素。v由由 有：有：得得 的第的第1行元素行元素 得得 的第的第1列元素到第列元素到第 列元素。由列元素。由 利用矩阵乘法有：利用矩阵乘法有：故故 又由又由 有：有：故：故：因此可得因此可得 的第的第 行和第行和第 列的全部元素。列的全部元素。直接分解法约需直接分解法约需 乘除法，和乘除法，和 消去法计算量基本消去法计算量基本相当。对计算相当。对计算 和式，可采用双精度累加以提高精度。和式，可采用双精度累加以提高精度。v（B）选主元的三角分解法）

34、选主元的三角分解法 在直接三角分解中，如果在直接三角分解中，如果 计算将要中断，或者当计算将要中断，或者当 绝对值很小时，按分解公式计算可能引起舍入误差的积累。但当绝对值很小时，按分解公式计算可能引起舍入误差的积累。但当 为非奇异时，可通过交换为非奇异时，可通过交换 的行实现矩阵的行实现矩阵 的的 分解。因分解。因此可采用与列主之消去法类似的方法将直接三角分解法修改为部分此可采用与列主之消去法类似的方法将直接三角分解法修改为部分选主之的三角分解法。选主之的三角分解法。设已完成第设已完成第 步分解，这时有：步分解，这时有：v第第 步分解需要到步分解需要到 两式，两式，为避免为避免 中中 用小的数

35、用小的数 作除数，先引进作除数，先引进量：量：由由 及及 定义，易见定义，易见 则由则由 有：有：若若 ，则我们可以用，则我们可以用 作为作为 交换交换 的第的第 行与行与 行（但我们将交换后的新元素仍行（但我们将交换后的新元素仍记为记为 及及 ）于是有）于是有 。控制了误差传播，再进行第。控制了误差传播，再进行第 步分解。步分解。对一般的非奇异矩阵对一般的非奇异矩阵 求逆的方法：求逆的方法：则：则：平方根法平方根法v 利用对称正定矩阵的三角分解而得到的求解利用对称正定矩阵的三角分解而得到的求解对称正定方程组的一种有效方法对称正定方程组的一种有效方法平方根法平方根法设设 为对称阵，即为对称阵，

36、即 ，且，且 的所有顺序的所有顺序主子式均非零，则知主子式均非零，则知 可以唯一分解为：可以唯一分解为：为利用为利用 的非奇异性，将的非奇异性，将 再分解为：再分解为：为对角矩阵，为对角矩阵，为单位上三角阵，则：为单位上三角阵，则：又又 ，由分解的唯一性即得：，由分解的唯一性即得：代入上面代入上面 中得：中得：。（对称阵的三角分解）设（对称阵的三角分解）设 为为 阶对称阵，且阶对称阵，且 的所有顺序主子式均非零。则的所有顺序主子式均非零。则 可唯一分解为：可唯一分解为：。其中。其中 为单位下三角阵，为单位下三角阵，为对角阵。为对角阵。当当 为对称正定矩阵时，则为对称正定矩阵时，则 的所有顺序主

37、子式的所有顺序主子式 ，而，而 设设 ，则，则 于是于是：从而：从而：其中其中 为下三角阵。为下三角阵。v定理定理 5.3.7 （对称正定矩阵的（对称正定矩阵的 分解）分解）若若 为为 阶对称正定矩阵，则存在一个实的非奇阶对称正定矩阵，则存在一个实的非奇异下三角矩阵异下三角矩阵 ，使，使 。当限定。当限定 的对角的对角元素为正时，这种分解是唯一的。元素为正时，这种分解是唯一的。下面来考虑计算下面来考虑计算 元素的公式，由元素的公式，由 由矩阵乘法及由矩阵乘法及 （时）时）有：有：。于是得解。于是得解对称正定方程组对称正定方程组 的平方根法计算公式：的平方根法计算公式：求解求解 即求解两个三角形

38、方程组即求解两个三角形方程组:(1)(2)由由 知知 ，则，则 ，从而，从而 。这表明分解过程中元素这表明分解过程中元素 的数量级不会增长太快且对的数量级不会增长太快且对角角元素元素 恒为整数。于是不选主元素的平方根是一个数值稳恒为整数。于是不选主元素的平方根是一个数值稳定的方法。定的方法。当求出当求出 的第的第 列元素时，列元素时，的第的第 行亦得出，所以行亦得出，所以平方平方根法大约需要根法大约需要 次乘除法，约为一般直接次乘除法，约为一般直接 分解法计分解法计算量算量的一半。由于的一半。由于 为对称阵，因此在计算机中只需存储为对称阵，因此在计算机中只需存储 的的下下三角部分元素，共需三角

39、部分元素，共需 个元素。可用一维数组存储，个元素。可用一维数组存储，即：即：。矩阵矩阵 的元素的元素 一维数组的表示为：一维数组的表示为：。的元素的元素存放在存放在 的相应位置。的相应位置。由公式由公式 可见，用平方根法解对称正定方程组时，可见，用平方根法解对称正定方程组时，计算计算 的元素的元素 需要进行开方计算，为避免开方计算，需要进行开方计算，为避免开方计算，可可对平方根法进行改进。对平方根法进行改进。用平方根法求解：用平方根法求解：精确解为精确解为 。(1)分解计算：分解计算：故，故，(2)求解两个三角方程：解求解两个三角方程：解 得得 代入代入 解得解得：。5.3.10 追赶法追赶法

40、v在一些实际问题中常有解三对角线性方程组在一些实际问题中常有解三对角线性方程组Ax=fAx=f的问题，的问题，即：即：其中其中 满足条件：满足条件：对于具有条件对于具有条件 的方程组的方程组 ，我们介绍下面的追，我们介绍下面的追赶法求解。赶法求解。追赶法追赶法具有计算量少，方法简单，算法稳定具有计算量少，方法简单，算法稳定的特点。的特点。v 设有三对角线性方程组设有三对角线性方程组Ax=fAx=f，且，且A A 满足条件满足条件 ，则，则A A为非奇异矩阵。为非奇异矩阵。v证明证明 用归纳法证明。显然用归纳法证明。显然n=2n=2时，有：时，有：v现假设定理对现假设定理对n-1n-1阶的满足条

41、件阶的满足条件 的三对角矩阵成立，的三对角矩阵成立，求证对满足条件求证对满足条件 的的 阶三对角矩阵定理亦成立。由阶三对角矩阵定理亦成立。由条件条件 ，则利用消元法的第，则利用消元法的第 步有：步有：显然显然 ，其中，其中 且有且有 故知故知 满足条满足条件件 ，利用归纳设知，利用归纳设知 ，故，故 设设 满足条件满足条件 为三对角阵。则为三对角阵。则 的所有的所有顺序主子式都不为零。即：顺序主子式都不为零。即：证明证明 由于由于A是满足条件（是满足条件（2）的）的n阶三对角阵。因此阶三对角阵。因此A的任的任一个顺序主子式亦满足条件（一个顺序主子式亦满足条件（2）的）的n阶三对角矩阵。由上阶三

42、对角矩阵。由上一个定理即知：一个定理即知：根据这一结论以及三角分根据这一结论以及三角分解定理知，这种矩阵解定理知，这种矩阵A可进行三角分解：可进行三角分解：A=LU。在这里特别的有：。在这里特别的有：其其中中待待定定系系数数 可可由由 利利用用矩阵乘法规则立即得出：矩阵乘法规则立即得出：由由 故故 v下面归纳证明：下面归纳证明：上面已经验证上面已经验证 对对 成立。现设成立。现设 对对 成立。求证成立。求证 对对 成立。成立。由假设由假设 ，再由，再由 及假设条件及假设条件 有有:故，故，由由 可知可知:可见由可见由A的假设条件的假设条件,我们完全求出了我们完全求出了 实现了对实现了对A的的L

43、U分解。从而求解分解。从而求解Ax=f 等价于求解两个三等价于求解两个三角方程组角方程组v由此可得求解三对角线性方程组的由此可得求解三对角线性方程组的追赶法追赶法:1.计算计算 的递推公式的递推公式:2.求解求解 Ly=f 的公式的公式:3.求解求解 Ux=y 的公式的公式:我们将计算我们将计算 及及 的过程称为追的过程。计算方程组的解的过程称为追的过程。计算方程组的解 的过程称为赶的过程。追赶法求解的过程称为赶的过程。追赶法求解 Ax=f 仅需仅需5n-4次乘次乘除运算。工作量较小。除运算。工作量较小。在计算机上，只需用在计算机上，只需用3个一维数组分别存储个一维数组分别存储A的系数的系数

44、以及两个一维数组保存计算的中间结果以及两个一维数组保存计算的中间结果 和和或或 。v例例5.3.7 用追赶法求解：用追赶法求解：v解解 计算计算 ：计算计算 ：计算计算 ：对于追赶法，由于已证明：对于追赶法，由于已证明：及及 因此追赶法中不会出现中间结果迅速增长和舍入误因此追赶法中不会出现中间结果迅速增长和舍入误差产生积累现象。差产生积累现象。5.4 小行星轨道方程问题的模型求解小行星轨道方程问题的模型求解v上上面面介介绍绍了了求求解解方方程程组组的的直直接接解解法法，针针对对此此题题，我我们们采采用三角分解法中一种方法求解，即用三角分解法中一种方法求解，即分解法分解法求解：求解：系数矩阵系数

45、矩阵 用用MATLABMATLAB软软件件求求解解，将将系系数数矩矩阵阵分分为为单单位位下下三三角角矩阵和上三角矩阵如下：矩阵和上三角矩阵如下：最后结果为：最后结果为：因此小行星的轨线方程为：因此小行星的轨线方程为：习习 题题 5 5v1 1、用高斯消去法求解如下线性方程组、用高斯消去法求解如下线性方程组v2 2、用主元素消去法求解如下方程组、用主元素消去法求解如下方程组v3 3、用列主元高斯消去法求解如下方程组、用列主元高斯消去法求解如下方程组v4 4、用直接三角分解法解如下方程组、用直接三角分解法解如下方程组v5 5、用、用DoolittleDoolittle分解法求方程分解法求方程 的解的解.v6 6、用平方根法求以下方程组的解、用平方根法求以下方程组的解.v 7 7、用法求解如下矩阵的逆矩阵、用法求解如下矩阵的逆矩阵v 8 8、用直接三角分解法求解如下方程组、用直接三角分解法求解如下方程组v 9 9、用追赶法求解如下方程组、用追赶法求解如下方程组中南大学数学科学与计算技术学院中南大学数学科学与计算技术学院