改进的平方根法及其程序实现.doc

1、改进的平方根法及程序实现目 录摘 要20 引 言31 预备知识31.1 分解理论31.2 Cholesky分解法 41.3 算法描述 5 2 改进的平方根法 63 分解算法描述 74 应用举例85 程序实现 105.1 程序码源 105.2 实例计算 126 结束语 13参考文献 14致 谢 15改进的平方根法及其程序实现摘要: 针对对称正定方程组的解法, 本文先对Cholesky分解法进行了分析研究,在此基础上给出了改进的平方根法(即分解法), 此方法有效地避免了原平方根法开方运算所带来的误差和不便, 并通过算法描述、实例计算, 用C程序实现了分解, 进一步提高了矩阵运算的速度和精度.关键词

2、: 对称正定矩阵, 平方根法, 分解, 算法Improved Methods of Square Root and Realization of Its Program Abstract: Aiming at studying solutions of symmetric positive definition matrix in linear equations. Initially, the text has conducted a series of analyses and researches towards decomposition proposed by Cholesky. T

3、hen based on theses researches and analyses, it offers the improved methods of squareroot (also called decom-position), which effectively avoid some errors and inconvenience brought by the process of extracting root. At the same time, it achieves the decomposition through the means of algorithm desc

4、ription, example calculation as well as applicat-ion of C program, further enhancing the speed and accuracy in matrix operation. Key words: Symmetric positive definition matrix, method of square root, decomposition, algorithm0 引言很多工程中的科学计算, 例如应用有限元法解结构力学问题时, 最后往往归结为求解系数矩阵为对称正定方程组解的问题. 由于对称正定矩阵各阶顺序主子

5、式以及全部的特征值均大于零, 这种特征也使得其三角分解具有更为简单的形式, 不同的分解也导出了一些不同的解法. 平方根法(即Cholesky分解法), 就是利用对称正定矩阵的三角分解而得到的求解对称正定方程组的一种有效方法, 其计算量和存储量约为普通消去法的一半, 且无需选主元就能求得较为精确的数值解, 但由于在平方根法中含有多次开方运算, 因此给计算带来了许多不便, 而在原平方根法的基础上, 给出改进的平方根法(即分解法), 成功避免了开方运算带来的的麻烦, 因此在各种工程计算中应用更为广泛.1 预备知识1.1 分解理论定理1.1（矩阵的分解定理） 设为阶矩阵, 如果的顺序主子式, 则可分解

6、为一个单位下三角矩阵和一个上三角矩阵的乘积, 且这种分解是唯一的.定理1.2（对称正定矩阵的分解定理） 设为阶对称正定矩阵, 则必存在非奇异下的三角矩阵, 使 , (1.1) 并且当的对角线元素均为正数时, 这种分解是唯一的.证明 因为对称正定, 则它的各阶顺序主子式均大于零. 由分解定理可知, 矩阵存在唯一的Doolitle分解, 即 其中, 为单位下三角矩阵, 为上三角矩阵. 记 diag, 则为单位上三角矩阵, 且. 由的对称性有, 可得 由于为单位下三角矩阵, , 从而由分解的唯一性可知, , 于是 显然有 ,这样 , 故可令 diag则 其中,为非奇异的下三角矩阵. 而知限定的对角元

7、素为正数时,这种分解是唯一的, 定理得证.1.2 Cholesky分解法设为对称正定矩阵, 由上面分解定理可知, 则存在一个实的非奇异下三角矩阵, 使, 当限定的对角线元素为正数时, 这时的分解是唯一的.设 (1.2)其中, , .由比较和的对应元素, 可求得的元素如下:由 , 得 , .假设的第列元素已经求得, 下面求的第列元素. 注意到 ,得计算公式, 对, 有 (1.3)由此可逐行求得的全部元素, 从而由及得到求解方程组的公式 (1.4)上述方法称为Cholesky分解法. 由于计算的对角线元素要作次开平方运算, 故Cholesky分解法又称平方根法.1.3 算法描述步1 输入对称正定矩

8、阵和右端向量; 步2 Cholesky分解: 对 计算: 步3 用向前消去法解下三角方程组: , 对 计算 ; 步4 解对角形方程组: 对 计算 ;步5 用回代法解上三角形方程组: , 对 计算 .不难验证Cholesky分解法的乘除计算总量约为, 为一般矩阵分解计算量的一半. 虽然如此, 但其增加的次开平方运算是非常不利的, 下面引出改进的平方根法分解法.2 改进的平方根法为了避免开方运算, 对作分解, 即, 则 / 提出矩阵的对角元素 / 由对称正定, 可得, 令 可证 即. (1.5)是对角元素为1的单位下三角矩阵.对矩阵作分解, 共计算个矩阵元素; 对称矩阵的分解, 只需计算个元素,

9、减少了近一半的工作量. 借助分解计算公式, 容易得到分解计算公式.设有分解形式其中.在分解中可套用分解公式, 只要计算下三角矩阵和的对角元素.计算中只需要保存的元素, 的行列的元素用的表示. 由于对称正定矩阵的各阶主子式均大于零, 直接调用分解公式可完成分解计算.对于,有 (1.6). (1.7)3 分解算法描述步1 输入方程组阶数、系数矩阵和常数项.步2 FOR TO ; ; FOR TO ; ;步3 解方程组的步骤从略.由, 解方程组可分为三步完成:(1)解方程组, 其中. (1.8)(2)解方程组, 其中. (1.9)(3)解方程组. . (1.10)可以看出, 改进后的分解乘除运算量约

10、为, 计算过程也无需开方运算, 使得其运算更加简单易行, 因此改进后的分解法相对于原平方根法具有更好的实效性和可行性.4 应用举例例 试用分解法求解方程组解 由, 可得的方程组, 令, 则.计算步骤:(1)将直接分解为, 求出;(2)求解方程组;(3)求解方程组.现有, , , , , , , 即 由可得解得 由有 解得 5 程序实现 5.1 程序码源应用分解法解线性方程组, 为了便于计算, 编写如下程序, 通过计算机来实现线性方程组的求解, 简洁方便. 具体程序如下: Purpose：分解法解方程组# include # include # define MAX_N 20 /(x_i,y_i

11、)的最大维数/int main()int n; int i,j,k; int mi;double mx,tmp; static double aMAX_NMAX_N,bMAX_N,xMAX_N,yMAX_N,zMAX_N; static double LMAX_NMAX_N,dMAX_N; /输入Ax=b的维数/ printf(n Input n value(dim of Ax=b):); scanf(%d,&n); if(nMAX_N) printf(The input n is larger then MAX_N,please redefine the MAX_N.n);return 1;

12、 if(n=0) printf(Please input a number between 1 and %d.n, MAX_N);return 1; /输入Ax=b的A矩阵/ printf(Now input the matrix a(i,j),i,j=0, %d:n, n-1);for(i=0;in;i+) for(j=0;jn;j+) scanf(%lf,&aij); /输入b矩阵/ printf(Now input the matrix b(i),i=0, %d:n, n-1); for(i=0;in;i+) scanf(%lf,&bi); for(i=0;in;i+) Lii=1; /

13、u矩阵对角元素为1/ for(k=0;kn;k+) dk=akk; /计算d_i/ for(j=0;j=k-1;j+) dk-=(Lkj *dj*Lkj); for(i=k+1;in;i+) /计算L矩阵/ Lik=aik; for(j=0;j=k-1;j+) Lik-=(Lij *dj*Lkj); Lik/=dk; for(i=0;in;i+) /求解Lz=b的z/ zi=bi; for(j=0;j=i-1;j+) zi-=Lij*zj;for(i=0;i=0;i-) /求解Ly=x/ xi=yi; for(j=i+1;jn;j+) xi-=Lji*xj; printf(Solve x_i=

14、n); /输出/ for(i=0;in;i+) printf(%fn,xi);return 0;5.2 实例计算 用分解法求解方程组:程序输入输出Input n value(dim of Ax=b): 3 Now input the matrix a(i,j),i,j=0,2: 16 4 8 4 5 -4 8 -4 22 Now input the matrix b(i),i=0,2: -4 3 10 Solve x_i= 2.000000 4.000000 -2.2500006 结束语从矩阵分解角度看, 分解法与消元法本质上没有多大的区别, 其基本思想还是通过等价变换将线性方程组化为结构简单

15、、易于求解的形式. 在实际问题中经常遇到对称正定矩阵方程组的求解问题, 对于这种具有特殊性质的系数矩阵, 采用改进的平方根法是一种行之有效的计算方法, 它是在普通消去法的基础上建立起来的, 算法思想虽有些复杂, 但其计算量约为普通消去法的一半, 且有效地避免了开方运算所带来的不便, 因此改进后的分解法相对于原平方根法具有更好的实效性和可行性. 并通过C程序求解使得这类问题变得更加简单易行, 在大量工程计算中有着广泛而重要的应用.参考文献1 刘元亮. Matlab在对称正定矩阵的改进平方根分解法中的应用J. 湖南: 怀化学院学报. 2004,(02).2 郭丽杰, 周硕, 秦万广. 对称矩阵的改

16、进Cholesky分解在特征值问题中的应用J. 吉林: 东北电力学院学报. 2003,(02).3 杜廷松, 沈艳君, 覃太贵. 数值分析及实验M. 北京: 科学出版社. 2007.4 肖筱南, 赵来军, 党林立. 数值计算方法与上机实习指导M. 北京: 北京大学出版社. 2004.5 马昌凤, 林伟川. 现代数值计算方法M. 北京: 科学出版社. 2008.6 张韵华, 奚梅成, 陈长松. 数值计算方法与算法M. 北京: 科学出版社. 2000.7 李庆扬, 王能超, 易大义. 数值分析(4版)M. 武汉: 华中科技大学出版社. 2006.8 严克明, 欧志英, 刘树群. 数值计算基础M.

17、甘肃: 甘肃人民出版社. 2006.9 Ferenc Szidsrovzky, Sidney Yokowitz. Principles and procedures of numerical analysisM. 上海: 上海科学技术文献出版社. 1982.10 Rainer Kress. Numerical analysisM. 北京: 世界图书出版社. 2003.致 谢本论文经过将近半年的努力即将告以尾声, 从论文的选题到最后完稿, 整个过程都经过了认真地考虑、仔细地查阅和细心地修改. 在此, 我首先要感谢我的导师郭晓斌老师, 不管是我的论文选题还是论文的撰写, 以及资料的查阅等各方面,

18、他都给了我莫大的帮助与启发, 尤其是在论文的几次修改过程中, 郭老师以他广博的专业学识、严谨的治学精神和他耐心的指导态度, 才使我的论文能顺利完成. 谨再次向郭老师致以崇高的敬意和诚挚的感谢. 也祝愿郭老师身体健康、工作顺利, 在学术研究上取得更加辉煌的成就, 为更多的学子导航.同时, 我也要对给予本论文参考文献的所有学术专家和老师致以真挚的谢意, 是他们出版的书籍与发表的学术论文给了我很大的启示与指导, 才将论文完成.其次, 在即将毕业之际, 我要感谢数信学院所有的老师对我的细心教育与培养, 让我在四年的学习生涯中不仅学到了扎实的专业知识, 而且他们的言传身教使我受益匪浅, 他们严谨的治学态度和耐心教导学生的博爱精神也是我永远学习的榜样, 并将积极影响着我今后的学习和工作.最后, 我还要再次感谢我的母校西北师范大学给了我发展的平台, 感谢各位老师四年来对我的耐心教导和栽培. 张登科 2011年5月 14