第4章线性方程组的数值解法优秀PPT.ppt

《第4章线性方程组的数值解法优秀PPT.ppt》由会员分享，可在线阅读，更多相关《第4章线性方程组的数值解法优秀PPT.ppt（119页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第4章线性方程组的数值解法现在学习的是第1页，共119页4.1引言引言在工程技术和科学研究中，很多科学计算的问题往往直接或间接地在工程技术和科学研究中，很多科学计算的问题往往直接或间接地归结为求解线性代数方程组。常见的线性代数方程组是方程个数和未归结为求解线性代数方程组。常见的线性代数方程组是方程个数和未知量个数相同的阶线性方程组，它的一般形式是知量个数相同的阶线性方程组，它的一般形式是其中，其中，、为常数，为常数，为待求的未知量为待求的未知量(4.1)现在学习的是第2页，共119页用矩阵形式表示是用矩阵形式表示是其中：其中：一般一般。当系数矩阵。当系数矩阵非奇异（行列式不为零），即非奇异（行

2、列式不为零），即时，方程时，方程组组(4.1)有惟一解。有惟一解。现在学习的是第3页，共119页对于线性代数方程组的解法大致可分为两类，即直接解法和对于线性代数方程组的解法大致可分为两类，即直接解法和迭代解法。直接解法是指在假设没有舍入误差的条件下，经过有迭代解法。直接解法是指在假设没有舍入误差的条件下，经过有限次算数运算就能求得方程组精确解的方法。限次算数运算就能求得方程组精确解的方法。然而实际计算中由于舍入误差的影响，这类方法也只能求得然而实际计算中由于舍入误差的影响，这类方法也只能求得近似解；迭代解法就是从一个已知的初始近似值开始，按一定的近似解；迭代解法就是从一个已知的初始近似值开始，

3、按一定的法则逐步求出解的各个更准确的近似值的方法，它是用某种极限法则逐步求出解的各个更准确的近似值的方法，它是用某种极限过程去逐步逼近精确解的方法。过程去逐步逼近精确解的方法。在本章中主要介绍高斯消去法、列主元高斯消去法、约当消在本章中主要介绍高斯消去法、列主元高斯消去法、约当消去法、三角分解法等直接解法以及雅可比、高斯去法、三角分解法等直接解法以及雅可比、高斯赛德尔和超松赛德尔和超松弛等迭代解法。弛等迭代解法。现在学习的是第4页，共119页4.2高斯（高斯（Gauss）消去法）消去法4.2.1 4.2.1 高斯消去法的基本思想高斯消去法的基本思想 高斯消去法是最古老的求解线性代高斯消去法是最

4、古老的求解线性代数方程组的方法之一，高斯消去法是消数方程组的方法之一，高斯消去法是消去法的一种特殊形式，它包括消元与回去法的一种特殊形式，它包括消元与回代两个过程。代两个过程。现在学习的是第5页，共119页设设 ，对行计算乘数，对行计算乘数用用 乘以第乘以第 个方程后加到第个方程后加到第 个方程到第个方程到第 个方程中，消去个方程中，消去 第个方程到第第个方程到第 个方程的未知数个方程的未知数 ，得到，得到现在学习的是第6页，共119页 只只要要设设 就就可可以以继继续续进进行行消消元元，直直到到经经过过 次次消消元元后后，将将线线性性方程组方程组(4.1)(4.1)化为化为(4.2)(4.2

5、)所示上三角方程组（以上计算过程称为消元过程）。所示上三角方程组（以上计算过程称为消元过程）。(4.2)(4.2)消消元元过过程程结结束束后后，只只要要设设 ，对对上上三三角角方方程程组组就就可可以以自自下下而而上上逐逐步回代，依次求得步回代，依次求得 ，即，即现在学习的是第7页，共119页4.2.2 4.2.2 实现高斯消去法的基本步骤实现高斯消去法的基本步骤（1 1）输入方程组的阶数输入方程组的阶数 n n ，系数矩阵，系数矩阵A A和右端常数矩阵和右端常数矩阵b b。（2 2）消元过程消元过程:设设 ，对，对 ，计算：，计算：（3 3）回代过程回代过程（4 4）输出方程组的解。输出方程组

6、的解。（5 5）结束。结束。现在学习的是第8页，共119页现在学习的是第9页，共119页例例4.1 4.1 用高斯消去法求解线性方程组，式中用高斯消去法求解线性方程组，式中现在学习的是第10页，共119页#define N 3/*N#define N 3/*N为方程组系数矩阵的阶数为方程组系数矩阵的阶数*/*/int Gauss(float aNN,float bN)int Gauss(float aNN,float bN)int i,j,k,flag=1;int i,j,k,flag=1;float t;float t;for(i=0;iN-1;i+)for(i=0;iN-1;i+)if(a

7、ii=0)if(aii=0)flag=0;break;flag=0;break;else else for(j=i+1;jN;j+)/*for(j=i+1;jN;j+)/*消元过程开始消元过程开始*/*/t=-aji/aii;t=-aji/aii;bj=bj+t*bi;bj=bj+t*bi;for(k=i;kN;k+)for(k=i;kN;k+)ajk=ajk+t*aik;ajk=ajk+t*aik;return(flag);return(flag);现在学习的是第11页，共119页void zg_matric(float aNN,float bN)/*void zg_matric(float

8、 aNN,float bN)/*输出增广矩阵输出增广矩阵*/*/int i,j;int i,j;for(i=0;iN;i+)for(i=0;iN;i+)for(j=0;jN;j+)for(j=0;j=0;i-)for(i=N-2;i=0;i-)xi=bi;xi=bi;for(j=i+1;jN;j+)for(j=i+1;jN;j+)xi=xi-aij*xj;xi=xi-aij*xj;xi=xi/aii;xi=xi/aii;for(i=0;iN;i+)/*for(i=0;iN;i+)/*输出方程组的解输出方程组的解*/*/printf(x%d=%11.7fn,i+1,xi);printf(x%d=

9、%11.7fn,i+1,xi);现在学习的是第13页，共119页程序运行结果：程序运行结果：0.100000 0.300000 0.500000 1.000000 0.100000 0.300000 0.500000 1.000000 3.000000 1.000000 -1.000000 5.000000 3.000000 1.000000 -1.000000 5.000000 3.000000 5.000000 4.000000 3.000000 3.000000 5.000000 4.000000 3.000000 x1=5.4583359 x1=5.4583359 x2=-6.5416

10、689 x2=-6.5416689 x3=4.8333344 x3=4.8333344现在学习的是第14页，共119页4.3 4.3 列主元高斯消去法列主元高斯消去法4.3.1 4.3.1 列主元高斯消去法的基本思想列主元高斯消去法的基本思想 高斯消去法，一直是在假设高斯消去法，一直是在假设 的情况下进行的，并且是按照方程组的情况下进行的，并且是按照方程组中各方程所给定的顺序进行的，所以这种方法又称为顺序消去法。但当中各方程所给定的顺序进行的，所以这种方法又称为顺序消去法。但当 时，消元过程就无法进行。另外，即使时，消元过程就无法进行。另外，即使 ，但其绝对值很小时，用它作，但其绝对值很小时，

11、用它作除数时，根据数值运算中除数时，根据数值运算中“用绝对值很小的数做除数，舍入误差会增大，而用绝对值很小的数做除数，舍入误差会增大，而且严重影响计算结果的精度且严重影响计算结果的精度”的原则，这种方法在一定程度上也具有局限性。为的原则，这种方法在一定程度上也具有局限性。为了克服这一局限性，产生了列主元高斯消去法。了克服这一局限性，产生了列主元高斯消去法。这种消去法在高斯消去法每次进行消元之前先选取绝对值最大的元素作为主元素这种消去法在高斯消去法每次进行消元之前先选取绝对值最大的元素作为主元素（位于主对角线上的在消去过程中用作除数的元素称为主元素，简称主元）进行消去。（位于主对角线上的在消去过

12、程中用作除数的元素称为主元素，简称主元）进行消去。现在学习的是第15页，共119页 列列主主元元高高斯斯消消去去法法，是是在在高高斯斯消消去去法法的的消消去去过过程程中中，第第 步步()增增加加选选主主元元操操作作，成成为为：首首先先在在第第 列列中中，从从 中中选选出出绝绝对对值值最最大大的的元元素素，经经过过行行交交换换(把把绝绝对对值值最最大大者者所所在在的的行行与与第第 行行交交换换)把把绝绝对对值值最最大大的的元元素素交交换换到到 的的单单元元中中；然然后后作作GaussGauss消消去去法法的的消消去去过过程程中中第第 步步，初初等等变变换换产产生生第第 列列的的 个个零零元元素素

13、。回回代代过过程程与与GaussGauss消去法的回代过程相同。消去法的回代过程相同。以上就是列主元高斯消去法的操作步骤。由于以上就是列主元高斯消去法的操作步骤。由于 ，可证，可证明明 ，中至少有一个元素不为零，从理论上保中至少有一个元素不为零，从理论上保证了列主元高斯消去法的可行性。证了列主元高斯消去法的可行性。现在学习的是第16页，共119页4.3.2 4.3.2 实现列主元高斯消去法的基本步骤实现列主元高斯消去法的基本步骤（1 1）输入方程组的阶数）输入方程组的阶数n n，系数矩阵，系数矩阵A A和右端常数矩阵和右端常数矩阵b b 。（2 2）列主元素）列主元素 对对 ，从，从 中选出绝

14、对值最大的元素中选出绝对值最大的元素 ，将，将 行和行和 行交换后，再作第行交换后，再作第 次消元操作。次消元操作。（3 3）消元过程：）消元过程：对对 ,计算：计算：现在学习的是第17页，共119页（4）回代过程：（5）输出方程组的解。（6）结束现在学习的是第18页，共119页现在学习的是第19页，共119页现在学习的是第20页，共119页例例 4.2 4.2 用列主元高斯消去法求解线性方程组用列主元高斯消去法求解线性方程组 ，式中，式中现在学习的是第21页，共119页void ColGauss(float aNN,float bN)void ColGauss(float aNN,float

15、 bN)float t,max_fab;float t,max_fab;int i,j,k,l;int i,j,k,l;for(i=0;iN-1;i+)for(i=0;iN-1;i+)l=i;l=i;max_fab=fabs(aii);max_fab=fabs(aii);for(j=i+1;jN;j+)/*for(j=i+1;jmax_fab)if(fabs(ajimax_fab)l=j;max_fab=fabs(aji);l=j;max_fab=fabs(aji);if(il)/*if(il)/*交换交换i i、l l两行两行*/*/t=bi;t=bi;bi=bl;bi=bl;bl=t;bl

16、=t;for(j=i;jN;j+)for(j=i;jN;j+)t=aij;aij=alj;alj=t;t=aij;aij=alj;alj=t;for(j=i+1;jN;j+)/*for(j=i+1;jN;j+)/*消元过程开始消元过程开始*/*/t=-aji/aii;t=-aji/aii;bj=bj+t*bi;bj=bj+t*bi;for(k=i;kN;k+)for(k=i;kN;k+)ajk=ajk+t*aik;ajk=ajk+t*aik;现在学习的是第22页，共119页void zg_matric(float aNN,float bN)/*void zg_matric(float aNN,

17、float bN)/*输出增广矩阵输出增广矩阵*/*/int i,j;int i,j;for(i=0;iN;i+)for(i=0;iN;i+)for(j=0;jN;j+)for(j=0;jN;j+)printf(%10f,aij);printf(%10f,aij);printf(%10f,bi);printf(%10f,bi);printf(n);printf(n);printf(n);printf(n);现在学习的是第23页，共119页#include#include#define N 4 /*N#define N 4 /*N为方程组系数矩阵的阶数为方程组系数矩阵的阶数*/*/main()m

18、ain()float aNN=1.003,0.333,1.504,-0.333,-2.011,1.455,0.506,2.956,float aNN=1.003,0.333,1.504,-0.333,-2.011,1.455,0.506,2.956,4.329,-1.952,0.006,2.087,4.329,-1.952,0.006,2.087,5.113,-4.004,3.332,-1.1125.113,-4.004,3.332,-1.112;float bN=3.005,5.407,0.136,3.772;float bN=3.005,5.407,0.136,3.772;float xN

19、=0,0,0,0;float xN=0,0,0,0;int i,j;int i,j;zg_matric(a,b);zg_matric(a,b);ColGauss(a,b);ColGauss(a,b);xN-1=bN-1/aN-1N-1;/*xN-1=bN-1/aN-1N-1;/*回代过程开始回代过程开始*/*/for(i=N-2;i=0;i-)for(i=N-2;i=0;i-)xi=bi;xi=bi;for(j=i+1;jN;j+)for(j=i+1;jN;j+)xi=xi-aij*xj;xi=xi-aij*xj;xi=xi/aii;xi=xi/aii;for(i=0;iN;i+)/*for(

20、i=0;iN;i+)/*输出方程组的解输出方程组的解*/*/printf(x%d=%11.7fn,i+1,xi);printf(x%d=%11.7fn,i+1,xi);现在学习的是第24页，共119页程序运行结果：程序运行结果：1.003000 0.333000 1.504000 -0.333000 3.005000 1.003000 0.333000 1.504000 -0.333000 3.005000 -2.011000 1.455000 0.506000 2.956000 5.407000 -2.011000 1.455000 0.506000 2.956000 5.407000 4.

21、329000 -1.952000 0.006000 2.087000 0.136000 4.329000 -1.952000 0.006000 2.087000 0.136000 5.113000 -4.004000 3.332000 -1.112000 3.772000 5.113000 -4.004000 3.332000 -1.112000 3.772000 x1=-0.325438 x1=-0.325438 x2=0.327990 x2=0.327990 x3=2.372718 x3=2.372718 x4=1.040163 x4=1.040163现在学习的是第25页，共119页 列主

22、元高斯消去法列主元高斯消去法,虽然保证了当前的主元素虽然保证了当前的主元素 是第是第 列列中中 以下元素中的绝对值最大者，但还不能保证它是同一行以下元素中的绝对值最大者，但还不能保证它是同一行（既第（既第 行）中的绝对值最大者。因此，经过列选主元后，虽然行）中的绝对值最大者。因此，经过列选主元后，虽然避免了避免了 的情况发生，但其计算过程还是不稳定的，不适于的情况发生，但其计算过程还是不稳定的，不适于求解较大规模的线性方程组，因而这种方法在一定程度上也具有求解较大规模的线性方程组，因而这种方法在一定程度上也具有局限性。为了克服这一局限性，产生了全主元高斯消去法。局限性。为了克服这一局限性，产生

23、了全主元高斯消去法。列主元高斯消去法当变换到第列主元高斯消去法当变换到第 步时，从系数矩阵的右下角步时，从系数矩阵的右下角 阶子矩阵中选取绝对值最大的元素，然后通过行变换与列变换将它交阶子矩阵中选取绝对值最大的元素，然后通过行变换与列变换将它交换到主元素换到主元素 的位置上。回代过程与的位置上。回代过程与GaussGauss消去法的回代过程相同，消去法的回代过程相同，请读者参考有关书籍。请读者参考有关书籍。现在学习的是第26页，共119页4.4 4.4 约当（约当（Jordan）消去法）消去法4.4.1 4.4.1 约当消去法的基本思想约当消去法的基本思想 约约当当（JordanJordan）

24、消消去去法法是是一一种种无无回回代代过过程程的的求求解解线线性性方方程程组组的的直直接接解解法法，基基本本思思想想是是，将将线线性性方方程程组组的的系系数数矩矩阵阵通通过过初初等等变变换换，化化为为一一个个等等价价的的对对角角阵阵（只只有有主主对对角角线线上上的的元元素不为素不为0 0，其余元素均为，其余元素均为0 0），从而方便地求出方程组的解。），从而方便地求出方程组的解。现在学习的是第27页，共119页 假如我们将此增广矩阵通过次初等变换，化成下列等价形式（式中每个元假如我们将此增广矩阵通过次初等变换，化成下列等价形式（式中每个元素右上角圆括号内的数字表示消去计算的次数）素右上角圆括号内

25、的数字表示消去计算的次数）由此可立即求得方程组的解：由此可立即求得方程组的解：现在学习的是第28页，共119页4.4.2 4.4.2 实现约当消去法的基本步骤实现约当消去法的基本步骤现在学习的是第29页，共119页现在学习的是第30页，共119页例例 4.3 4.3 用约当消去法求解线性方程组用约当消去法求解线性方程组 ，式中，式中，现在学习的是第31页，共119页#define N 3 /*N#define N 3 /*N为方程组的系数矩阵的阶数为方程组的系数矩阵的阶数*/*/int joydan(float aNN+1)int joydan(float aNN+1)float r;floa

26、t r;int i,j,k,flag=1;int i,j,k,flag=1;for(k=0;kN;k+)for(k=0;kN;k+)if(akk=0)/*if(akk=0)/*当当akkakk为为0 0时，时，约当消去法不能进行约当消去法不能进行*/*/flag=0;flag=0;break;break;else else for(i=0;iN;i+)for(i=0;iN;i+)r=aik/akk;r=aik/akk;if(i!=k)if(i!=k)if(aik!=0)if(aik!=0)for(j=0;j=N;j+)aij=aij-r*akj;for(j=0;j=N;j+)aij=aij-r

27、*akj;return(flag);return(flag);现在学习的是第32页，共119页void print_matric(float aNN+1)/*void print_matric(float aNN+1)/*输出增广矩阵输出增广矩阵*/*/int i,j;int i,j;for(i=0;iN;i+)for(i=0;iN;i+)for(j=0;j=N;j+)for(j=0;j=N;j+)printf(%12f,aij);printf(%12f,aij);printf(n);printf(n);printf(n);printf(n);现在学习的是第33页，共119页main()mai

28、n()float aNN+1=1,-1,1,-4,5,-4,3,-12,2,1,1,11;float aNN+1=1,-1,1,-4,5,-4,3,-12,2,1,1,11;float xN;float xN;int i,k;int i,k;clrscr();clrscr();print_matric(a);/*print_matric(a);/*输出增广矩阵的初始数据输出增广矩阵的初始数据*/*/k=joydan(a);k=joydan(a);if(k=1)/*if(k=1)/*当当k=1k=1时时约当消去法正常进行约当消去法正常进行*/*/print_matric(a);/*print_m

29、atric(a);/*输出初等变换后的增广矩阵输出初等变换后的增广矩阵*/*/for(i=0;iN;i+)for(i=0;iN;i+)xi=aiN/aii;xi=aiN/aii;printf(x%d=%11.6fn,i+1,xi);/*printf(x%d=%11.6fn,i+1,xi);/*输出方程组的解输出方程组的解*/*/else /*else /*当当k=0k=0时输出时输出约当消去法不能进行的信息约当消去法不能进行的信息*/*/printf(Joydan method does not run!n);printf(Joydan method does not run!n);现在学习的

30、是第34页，共119页程序运行结果：程序运行结果：1.000000 -1.000000 1.000000 -4.000000 1.000000 -1.000000 1.000000 -4.000000 5.000000 -4.000000 3.000000 -12.000000 5.000000 -4.000000 3.000000 -12.000000 2.000000 1.000000 1.000000 11.000000 2.000000 1.000000 1.000000 11.000000 1.000000 0.000000 0.000000 3.000000 1.000000 0.

31、000000 0.000000 3.000000 0.000000 1.000000 0.000000 6.000000 0.000000 1.000000 0.000000 6.000000 0.000000 0.000000 5.000000 -5.000000 0.000000 0.000000 5.000000 -5.000000 x1=3.000000 x1=3.000000 x2=6.000000 x2=6.000000 x3=-1.000000 x3=-1.000000 约当消取法能否进行，取决于每一步能否满足约当消取法能否进行，取决于每一步能否满足 ，否则约当消取法将无，否则约

32、当消取法将无法进行。另外，即使法进行。另外，即使 ，但其绝对值很小时，仍用它作除数时，就会导致，但其绝对值很小时，仍用它作除数时，就会导致舍入误差增大，严重影响计算结果的精度。为解决这一问题，可以采用选主元的技舍入误差增大，严重影响计算结果的精度。为解决这一问题，可以采用选主元的技术进行计算，具体算法可参考有关书籍。术进行计算，具体算法可参考有关书籍。现在学习的是第35页，共119页4.5 4.5 三角分解法三角分解法 三三角角分分解解法法是是高高斯斯消消去去法法的的一一种种变变形形解解法法。用用高高斯斯消消去去法法去去解解 阶阶线线性性方方程程组组 ，经经过过 次次消消元元之之后后，得得出出

33、一一个个等等价价的的上上三三角角形形方方程程组组 ，对对上上三三角角形形方方程程组组用用逐逐步步回代的方法就可以求出解来。回代的方法就可以求出解来。上述消元过程也可以通过矩阵分解来实现。上述消元过程也可以通过矩阵分解来实现。现在学习的是第36页，共119页4.5.1 4.5.1 三角分解法的基本思想三角分解法的基本思想 设设 阶阶线线性性方方程程组组(4.1)(4.1)的的系系数数矩矩阵阵 非非奇奇异异，并并且且系系数数矩矩阵阵 ，能分解为两个三角形矩阵之积，即，能分解为两个三角形矩阵之积，即 (4.3)(4.3)其其中中 为为一一个个下下三三角角阵阵，为为一一个个上上三三角角阵阵，那那么么方

34、方程程组组 就化为就化为 (4.4)(4.4)从而使问题转化为求解两个简单的三角形方程组：从而使问题转化为求解两个简单的三角形方程组：(4.5)(4.5)现在学习的是第37页，共119页 矩矩阵阵 分分解解为为两两个个矩矩阵阵 、之之积积，一一般般不不是是惟惟一一的的，但但规规定定了了 、的的形形式式之之后后，由由矩矩阵阵的的乘乘法法运运算算法法则则知知，这这种种分分解解便便是是惟惟一一的的。由由于于 和和 、的的形形式式不不同同，就就得得到到矩矩阵阵的不同解法，相应地就有解线性代数方程组的几种解法的不同解法，相应地就有解线性代数方程组的几种解法.本本节节介介绍绍直直接接三三角角分分解解法法、

35、乔乔里里斯斯基基（CholeskyCholesky）分分解解法、追赶法和对称高斯法（法、追赶法和对称高斯法（分解法）。分解法）。现在学习的是第38页，共119页4.5.2 4.5.2 直接三角分解法直接三角分解法 设线性方程组设线性方程组(4.1)(4.1)的系数矩阵非奇异，并且系数矩阵，的系数矩阵非奇异，并且系数矩阵，能分解为同阶的单位下三角形矩阵和同阶的上三角矩阵之能分解为同阶的单位下三角形矩阵和同阶的上三角矩阵之积，或能分解为同阶的下三角形矩阵和同阶的单位上三角积，或能分解为同阶的下三角形矩阵和同阶的单位上三角矩阵之积，即写成矩阵元素的形式如下：矩阵之积，即写成矩阵元素的形式如下：现在学

36、习的是第39页，共119页或 则称为直接三角分解法。将矩阵 分解为 时，当 是单位下三角形矩阵，是上三角矩阵时，称为Dolittle分解；当 是下三角形矩阵，是单位上三角矩阵时，称为Courant分解。Dolittle分解法和Courant分解法的基本思想相同，本章只介绍Courant分解法。现在学习的是第40页，共119页4.5.3 4.5.3 实现实现CourantCourant分解法的基本步骤分解法的基本步骤现在学习的是第41页，共119页现在学习的是第42页，共119页现在学习的是第43页，共119页#define N 3 /*N#define N 3 /*N为方程组系数矩阵的阶数为方

37、程组系数矩阵的阶数*/*/void Lu(float aNN,float LNN,float UNN)/*LUvoid Lu(float aNN,float LNN,float UNN)/*LU分解分解*/*/int i,j,k;int i,j,k;for(k=0;kN;k+)for(k=0;kN;k+)for(i=k;iN;i+)/*for(i=k;iN;i+)/*计算计算L L矩阵的第矩阵的第k k列元素列元素*/*/Lik=aik;Lik=aik;for(j=0;j=k-1;j+)for(j=0;j=k-1;j+)Lik-=(Lij*Ujk);Lik-=(Lij*Ujk);for(j=k

38、+1;jN;j+)/*for(j=k+1;jN;j+)/*计算计算U U矩阵的第矩阵的第k k列元素列元素*/*/Ukj=akj;Ukj=akj;for(i=0;i=k-1;i+)for(i=0;i=k-1;i+)Ukj-=(Lki*Uij);Ukj-=(Lki*Uij);Ukj/=Lkk;Ukj/=Lkk;现在学习的是第44页，共119页void solve(float bN,float LNN,float UNN,float xN)void solve(float bN,float LNN,float UNN,float xN)int j,k;int j,k;float yN;float

39、yN;for(k=0;kN;k+)/*for(k=0;kN;k+)/*计算计算Ly=bLy=b中的中的y*/y*/yk=bk;yk=bk;for(j=0;j=k-1;j+)for(j=0;j=0;k-)/*for(k=N-1;k=0;k-)/*计算计算Ux=yUx=y中的中的x*/x*/xk=yk;xk=yk;for(j=k+1;jN;j+)for(j=k+1;jN;j+)xk=xk-Ukj*xj;xk=xk-Ukj*xj;现在学习的是第45页，共119页void zg_matric(float aNN,float bN)/*void zg_matric(float aNN,float bN)

40、/*输出增广矩阵输出增广矩阵*/*/int i,j;int i,j;for(i=0;iN;i+)for(i=0;iN;i+)for(j=0;jN;j+)for(j=0;jN;j+)printf(%10f,aij);printf(%10f,aij);printf(%10f,bi);printf(%10f,bi);printf(n);printf(n);printf(n);printf(n);现在学习的是第46页，共119页main()main()static float aNN=1,2,-1,1,-1,5,4,1,-2;static float aNN=1,2,-1,1,-1,5,4,1,-2;

41、float bN=3,0,2;float bN=3,0,2;float xN=0,0,0;float xN=0,0,0;float UNN=1,0,0,0,1,0,0,0,1;float UNN=1,0,0,0,1,0,0,0,1;float LNN;float LNN;int i;int i;zg_matric(a,b);zg_matric(a,b);Lu(a,L,U);Lu(a,L,U);solve(b,L,U,x);solve(b,L,U,x);for(i=0;iN;i+)/*for(i=0;iN;i+)/*输出方程组的解输出方程组的解*/*/printf(x%d=%11.7fn,i+1

42、,xi);printf(x%d=%11.7fn,i+1,xi);现在学习的是第47页，共119页程序运行结果：程序运行结果：1.000000 2.000000 -1.000000 3.000000 1.000000 2.000000 -1.000000 3.000000 1.000000 -1.000000 5.000000 0.000000 1.000000 -1.000000 5.000000 0.000000 4.000000 1.000000 -2.000000 2.000000 4.000000 1.000000 -2.000000 2.000000 x1=0.2500000 x1=

43、0.2500000 x2=1.5000000 x2=1.5000000 x3=0.2500000 x3=0.2500000现在学习的是第48页，共119页 在在直直接接三三角角分分解解法法中中计计算算时时，如如果果分分母母（甚甚至至可可能能为为零零）绝绝对对值值很很小小，仍仍用用它它作作除除数数时时，就就会会导导致致舍舍入入误误差差增增大大，严严重重影影响响计计算算结结果果的的精精度度。为为解解决决这这一一问问题题，可可以以配配合合按按列列选选主主元元的的技技术术进进行行计计算算，具具体体算算法法可可参参考考有有关关书书籍。籍。现在学习的是第49页，共119页4.5.4 4.5.4 乔里斯基（

44、乔里斯基（CholeskyCholesky）分解法）分解法 在在工工程程和和科科学学计计算算中中，有有时时归归结结为为求求一一个个特特殊殊的的线线性性代代数数方方程程组组，例例如如系系数数矩矩阵阵是是一一个个对对称称正正定定矩矩阵阵，为为此此，这这里里介介绍绍系系数数矩矩阵阵为为对对称称正正定定矩矩阵阵的的方方程程组组的的三三角角解解法法乔乔里里斯斯基基分分解法。解法。现在学习的是第50页，共119页 对对 阶阶线线性性方方程程组组 ，当当系系数数矩矩阵阵 为为对对称称正正定定时时可可以以惟惟一一地地分分解解为为 ，其其中中 为为上上三三角矩阵。即写成矩阵元素的形式如下：角矩阵。即写成矩阵元素

45、的形式如下：且且 。现在学习的是第51页，共119页4.5.5 4.5.5 实现乔里斯基分解法的基本步骤实现乔里斯基分解法的基本步骤现在学习的是第52页，共119页现在学习的是第53页，共119页实实现现乔乔里里斯斯基基分分解解法法的的N NS S图图现在学习的是第54页，共119页现在学习的是第55页，共119页例例 4.5 4.5 用乔里斯基分解法求解线性方程组，式中用乔里斯基分解法求解线性方程组，式中现在学习的是第56页，共119页#include#include#define N 4 /*N#define N 4 /*N为方程组系数矩阵的阶数为方程组系数矩阵的阶数*/*/void ch

46、odec(float aNN,float lNN)void chodec(float aNN,float lNN)int i,j,k;int i,j,k;float s;float s;for(i=0;iN;i+)li0=ai0/sqrt(a00);for(i=0;iN;i+)li0=ai0/sqrt(a00);for(i=1;iN;i+)for(i=1;iN;i+)s=0.0;s=0.0;for(k=0;k=i-1;k+)for(k=0;k=i-1;k+)s=s+lik*lik;s=s+lik*lik;lii=sqrt(aii-s);lii=sqrt(aii-s);for(j=i+1;jN;

47、j+)for(j=i+1;jN;j+)s=0.0;s=0.0;for(k=0;k=i-1;k+)for(k=0;k=i-1;k+)s=s+ljk*lik;s=s+ljk*lik;lji=(aji-s)/lii;lji=(aji-s)/lii;现在学习的是第57页，共119页void solve(float aNN,float bN,float xN)void solve(float aNN,float bN,float xN)int i,j,k;int i,j,k;float s;float s;for(k=0;kN-1;k+)for(k=0;kN-1;k+)bk/=akk;bk/=akk;f

48、or(i=k+1;iN;i+)for(i=k+1;i=0;k-)for(k=N-2;k=0;k-)s=bk;s=bk;for(j=k+1;jN;j+)for(j=k+1;jN;j+)s-=ajk*xj;s-=ajk*xj;xk=s/akk;xk=s/akk;现在学习的是第58页，共119页void zg_matric(float aNN,float bN)void zg_matric(float aNN,float bN)/*/*输出增广矩阵输出增广矩阵*/*/int i,j;int i,j;for(i=0;iN;i+)for(i=0;iN;i+)for(j=0;jN;j+)for(j=0;j

49、N;j+)printf(%10f,aij);printf(%10f,aij);printf(%10f,bi);printf(%10f,bi);printf(n);printf(n);printf(n);printf(n);现在学习的是第59页，共119页main()main()static float aNN=1,2,1,-3,2,5,0,-5,1,0,14,1,static float aNN=1,2,1,-3,2,5,0,-5,1,0,14,1,-3,-5,1,15;-3,-5,1,15;float lNN;float lNN;float bN=1,2,16,8;float bN=1,2,

50、16,8;float xN=0,0,0,0;float xN=0,0,0,0;int i;int i;zg_matric(a,b);zg_matric(a,b);Chodec(a,l);/*Chodec(a,l);/*分解过程分解过程*/*/Solve(l,b,x);/*Solve(l,b,x);/*回代求方程组的解回代求方程组的解*/*/for(i=0;iN;i+)/*for(i=0;iN;i+)/*输出方程组的解输出方程组的解*/*/printf(x%d=%11.7fn,i+1,xi);printf(x%d=%11.7fn,i+1,xi);现在学习的是第60页，共119页程序运行结果：程序