蒙伟拟Toeplitz带状矩阵线性方程组的并行算法(最终版)

1、拟Toeplitz带状矩阵线性方程组的并行算法蒙 伟（咸阳师范学院 数学与信息科学学院 陕西 咸阳 712000 ）摘 要本文主要研究了Toeplitz矩阵线性方程组的并行算法。首先介绍了Toeplitz矩阵的定义、分类及半正定性等，然后分别给出了三对角Toeplitz线性方程组、拟对称七对角Toeplitz线性方程组的并行算法的误差分析。最后，用Matlab语言编写了算法程序，通过实例计算，结果表明本文算法是可行的，并且验证了理论分析与实际计算的结果是一致的。关键词：Toeplitz矩阵；Toeplitz线性方程组；误差分析；并行算法；并行性分析；拟对称七对角Toeplitz线性方程组Par

2、allel Algorithm for Solving Banded quasi-Toeplitz Matrix Linear EquationsMeng wei(School of Mathematics and Information Science of Xianyang NormalUniversity Xianyang 712000)AbstractThis paper studies parallel algorithm of the Toeplitz matrix linear equations. Firstly, I introduce the definition, cla

3、ssification, and semi-positive definite of Toeplitz matrix and so on, and then give the error analysis of parallel algorithm for solving the three-diagonal Toeplitz matrix linear equations and quasi-symmetry seven diagonal Toeplitz matrix linear equations. Finally, I program the algorithm procedure

4、using Matlab software. The results show that the algorithm is feasible, and the theoretical analysis is consistent with the result of practical calculation. KEY WORDS：Toeplitz matrix; Toeplitz matrix linear equations; Error analysis; Parellel algorithm; Quasi-symmetry seven diagonal Toeplitz matrix

5、linear equations目录摘要IABSTRACTII前言11 TOEPLITZ矩阵及其性质31.1 Toeplitz的定义31.2 Toeplitz的半正定性42 三对角TOEPLITZ矩阵线性方程组的并行算法52.1 并行算法52.2 误差分析82.3 并行性分析92.4 算例103 拟对称七对角TOEPLITZ矩阵线性方程组的并行算法113.1并行算法113.2 误差分析164 结论17参考文献18附录19谢辞24前 言Toeplitz矩阵是数学和工程中应用广泛的特殊矩阵之一。由于工程中许多典型问题都涉及Toeplitz线性方程组的求解，所以有必要就Toeplitz线性方程组的求

6、解进行详细介绍。利用Toeplitz矩阵具有的非正常特殊的结构，我们可以由三对角Toeplitz矩阵线性方程组的并行算法得出对称七对角Toeplitz矩阵线性方程组的并行算法。随着高性能并行机的出现，为了更好的利用并行机，研究并行算法已知变的越来越重要了。在许多的工程问题中出现Toeplitz矩阵的线性方程组的求解问题，研究该问题的并行算法具有重要的现实意义。目前，对于对称三对角Toeplitz线性方程组的并行算法的研究，已获得了良好的结果，并利用此推广到了对称五对角Toeplitz线性方程组的求解，在此我先对L.EGarey的研究过程作简要的介绍。首先，我们先看看L.EGarey对五对角的处

7、理过程由 ，可得，计算 进行LU分解得：，且在中，是一个阶的形式的矩阵，是阶零矩阵。在中，是一个阶的形式的矩阵，是只有其它位置为0的矩阵。只有4个非零元素：，即 （a） （b）误差分析：这就是L.EGarey对五对角线的主要处理过程。在本文中，我们利用此思想，在文章的第一章中给出Toeplitz的定义和其性质，其次在第二章介绍三对角Toeplitz矩阵的线性方程组的算法，在第三章给出特殊对称七对角线性方程组的算法，并给出了误差分析。通过上机计算表明，理论与实际是一致的。1 Toeplitz矩阵及其性质1.1 Toeplitz的定义在数学和工程（例如信号处理与系统理论等）问题中，常常需要求解具有

8、特殊结构的线性方程组，其中 （1.1）这种形式取的矩阵称为Toeplitz矩阵。因此，任何一条对角线取相同元素的矩阵是Toeplitz矩阵，之所以称为Toeplitz矩阵，乃是这种特殊矩阵是Toeplitz于1900年代初在研究与Laurent级数有关的双线性函数的一篇论文中提出来的。最常见的Toeplitz矩阵为对称Toeplitz矩阵，即其元素还满足对称关系，可见，对称Toeplitz矩阵仅由其第1行元素就可完全描述，因此，常将对称Toeplitz矩阵简记作。若一个复Toeplitz矩阵的元素满足复共轭对称关系，即 （1.2）则称之为Hermitian Toeplitz矩阵特别地，具有特殊

9、结构 （1.3）的维Toeplitz矩阵称为Hermitian Toeplitz矩阵；而 （1.4）称为Hermitian Toeplitz矩阵。式中，为实数，显然，一个斜Hermitian型Toeplitz矩阵可以写作，其中，为斜Hermitian Toeplitz矩阵，斜Hermitian Toeplitz矩阵和斜Hermitian 型Toeplitz矩阵在求解离散化的双曲差分方程时经常会遇到。1.2 Toeplitz的半正定性Toeplitz矩阵在数学和工程中有着广泛的应用。例如，在信号处理中，有赖于通过求解Toeplitz矩阵方程组获得的参数包括：递推数字滤波器的系数，一维和二维平稳回

10、归（AR）模型的AR参数等，又如，信号和图像的恢复也可归结为Toeplitz矩阵方程组的求解，其他的应用例子有：偏微分方程和卷积型方程的求解，Pade逼近和控制理论中的最小实现问题等。考查式（1.1）的对称Toeplitz矩阵，令其主子式为， （1.5）则是正定的，当且仅当，然而，值得指出的是，类似的结论对于对称Toeplitz矩阵的半正定性却不一定成立。更具体地说，这一事实对于的半正定性并不是充分条件，请看一个反例： （1.6）在这个例子中，但是，删去矩阵的第2行和第2列组成的主子式都不大于或等于零时，对称Toeplitz矩阵才是半正定的，显然，这一条件的判断比较麻烦，下面给出了对称Toep

11、litz矩阵半正定性的一种简单检验方法，它不需要计算任何主子式。令，是一个对称Toeplitz矩阵，若是满足正定和条件的最小正整数，则矩阵，是半正定的，当且仅当系数服从递归方程， （1.7）式中，为模型的系数。Toeplitz矩阵具有以下性质：（1）Toeplitz矩阵的线性组合仍然为Toeplitz矩阵；（2）若Toeplitz矩阵的元素，则为对称Toeplitz矩阵；（3）Toeplitz矩阵的转置仍然为Toeplitz矩阵；（4）Toeplitz矩阵的元素相对于交叉对角线对称。2 三对角Toeplitz矩阵线性方程组的并行算法2.1 并行算法设方程组 （2.1）其他当时，（1）转化成 （

12、2.2）设，将进行如下分裂： （2.3）其中 ，其中，且，。可分解为：其中 （2.4）由（2.4）知取的解设其中 ，即：首先求解方程组：，显然可以化为两个方程组：这就把一个大方程组分解成两个小的方程组求解。如此类推我们可以把一个大方程组分解成3个，4个甚至更多，从而方便了方程组的并行求解。由，我们可以推得所以只需求解以下三个方程：，按照文献1中的方法取其中为中的根，即。由于所以因而解的近似解 （2.5）通过如上的推理我们得到了近似解的表达式（2.5）。2.2 误差分析在此，我们进行如下的误差分析：因为取向量的无穷范数，那么有我们可以得到如下结论。定理：由（2.4）式得到的近似解，那么满足由此可

13、见，得好坏直接取决于的大小，若要求，只需，由此可得到： （2.6）由于，所以若要达到精度要求，矩阵的阶数应大一些。当矩阵中时就是一个对称的。只要把用进行替代就可以处理三对角对称情形。2.3 并行性分析为了方便起见，我们不妨设，是处理机台数，每台处理机求解需要的时间为： 在修正过程需要时间运算量为： 并行通讯二次需要时间为： 并行计算时间为：串行元素按时间为： 所以 并行加速比为： 并行效率 按照并行算法的意义在时，可见该并行算法具有良好的并行性。2.4 算例方程中，为一120阶非对称Toeplitz矩阵我们求得方程组的近似解要求，由（2.6）式可得取实际计算结果得：我们举第二个例子如下：要求，

14、由（2.6）式可得取实际计算解果得：近似解和精确解的误差在分裂后生成的误差都小于这就意味着解的各个分量误差均小于，解的符合性较好。其中为把矩阵进行LU分解时的下三角矩阵中的元素，为非对称矩阵修正时的，是修正方程的近似解。由程序的运行结果和理论计算非常的接近，满足误差分析的精度要求。由于实际的计算结果和理论结果结合的比较合适，没有出现较大的误差。由于矩阵的元素取值非常的苛刻，对于随机输入元素的计算结果离理论结果非常的远，如主对角元素为，其出现的误差达到了，离我们的理论结果相去甚远。3 拟对称七对角Toeplitz矩阵线性方程组的并行算法3.1并行算法在许多的工程问题中都会出现特征值的求解问题。在

15、对特征值求解的研究中，向后差分法和有限分法都存在比较狭窄的使用前提。对于满足有边界条件和辅助条件的系统是完全的且系数矩阵是对称的，而以自然辅助条件取代边界条件和辅助条件的系统，带状系数矩阵变成了非对称的。对于带状拟对称系统引出了新的研究工作。为了寻找一个有效的并行解法求解工程中更为一般的问题，在文献1中引入了一个五角对角线拟Toeplitz带状对称矩阵的有效并行算法。为推广该方法的应用，我认为可以对其作进一步的推广，以适应更为复杂的问题的求解。本文以七对角线为例，对于更高阶的推广可以依次进行类型的拓展。首先，我们引入一个阶的系统 （3.1）其中 （3.2）矩阵就是一个对称七对角拟Toeplit

16、z矩阵。令 有令，则可写成如下形式参照文献1中五对角线矩阵的处理方法，将分解为其中为三对角线对称矩阵，为无对角线对称矩阵，分别为：由知道，我们可以得到如下等式有方程组 我们求得的值。由的值我们可以求得的值。在该处理过程中，需要满足为主对角占优矩阵，则必须成立。为一无对角线带状对称矩阵，则对可进行如同的分解，则，和均为一个三对角线带状对称矩阵。此时令此时或者其中必须有一个成立。不失一般性，不妨令，综合对的处理条件，我们要进行操作的充分条件是或者或者由于为对角占优（也是一样处理）可以写成如下形式其中其中，按照在文献2中的方法对进行分解得，我们注意对和都是对角占优的。由和的两个关系，令且 由这些值，

17、分解到和只需要少量的计算量。我们也可以把写出如下形式 （3.3），若，则，若，则在中，是一个阶的形式的矩阵，是阶零矩阵。在中，是一个阶的形式的矩阵，是只有其它位置为0的矩阵。只有4个非零元素：，即 右边的数字表示对应的行数由，得即：用文献2中的并行方法进行求解的近似解 （3.4）再求解，得（7）式的近似值，以此类推，这种方法可以广泛的推广到把系统分成4个或更多的子系统而用并行算法求解。除了一些明显的限制外，修正部分要求每个方程组的阶数足够的大，以适应最终近似的准确度要求。利用1中的方法，得到方程组的近似解，即就是（3.2）的近似解。要求（3.1）的解，必须对解进行修正。因为所以其中这时我们需要

18、求两个向量，使得，设，其中是方程个不同的根。显然，的确定只需的前三与后三个方程解得。而，这样我们得到了方程组（3.1）的近似解3.2 误差分析设， 分别是，的近似解，因为 =所以可见其误差直接取决于三个三对角线性方程组的解的误差，其近似度等于 所以解七对角线性方程组结果的好坏，关键的问题就是其近似分解出的每个三对角线性方程组的结果的好坏，只要三对角线性方程组的结果好，则七对角线性方程组的结果也不会差。我们只给出算法分析，由于方程最终的近似需要修正很多的误差，在程序上进行实现比较困难，我们用并行算法求解七对角甚至更高的时候，修正的部分比较多，为了实现更加近似的近似解，我们要对误差的修正做更加精确

19、的修正。4 结论本文采用并行数值解法解决了具有Toeplitz系数矩阵的七对角线性方程组，我们在论文的开头引入了我们要解决的系统并对其进行了第一步的修正处理，在论文的第二部分，我们具体的分析了对称系统的性质及其可以修正的修正过程，我们分别对第一部分中引入的系统作了LU分解和并行算法的处理，并比较了他们处理的条件和最后结算量在不同条件下的比较。在本文的第三部分，我们对扰动系统作了修正，并且找到第一部分引入系统的近似解，通过引入的定理对系统的最终解做出了近似并分析了各自在给定精度下的允许误差限。并最终求解出系统的最终近似解。在论文的最后一部分，我们再一次对处理系统的LU方法和本文提供的近似并行解法

20、的运算量作了最后的比较，并指出了近似并行解法的优越性和前景分析。参 考 文 献1 LEGarey A parallel numerical algorithm for near symmetric and banded systems 2 WMYanKLChung，A fast algorithm for soving special tridiagonal systems，Computing 52(1994)203-2113 JMMcNally，LEGarey，REShaw，A split-correct parallel algorithm for solving tridiagonal

21、symmetric Toeplitz systems，IntJComputMaths 4M.M.chawla，C.P.Katti，A new fourth-orded method for computing eigenvalue，of two-point boundary value problems，BIT20(1980)511-5145M.Chen，On the solution of circulant liner systems，SIAMJ.Numer.Anal.24(1987)668-6836L.E.Garey，R.E.Shaw，A parallel algorith for so

22、lving Toeplitz liner systems，Appl.Math.Comp.100(1999)241-2477J.H.Bramble，B.E.Hubbard，On afinite difference analogue of an elliptic boundary value problem which is nerthei diagonally dominant nor of non-negative type，J.Math.Phys.43(1964)117-1328A.Ralston，A First Course in Numerical Analysis，McGraw-Hi

23、ll，New York，19659O.Rojo，A new method for solving symmetrical circulant tridagonal systems of liner equations，Comput.Math.Appl.20(1990)61-6710 R.E.Shaw，L.E.Garey，A fast method for solving second order boundary value Volterra integro-differential equations，Int.J.Comput.Math 65(1997)121-12911 成田诚之助著.张贤

24、达译.数学系统控制理论及应用.北京：机械工业出版社，198412 陈景良，陈向辉.特殊矩阵.北京：清华大学出版社，200113 黄琳.系统与控制理论中的线性代数.北京：科学出版社，198414 数学手册编写组，数学手册.北京：高等教育出版社，197915 张贤达.现代信号处理.北京：清华大学出版社，199516 张贤达.现代信号处理（第二版）.北京：清华大学出版社，200217 张贤达.信号处理中的线性代数.北京：科学出版社，199718 张贤达.矩阵分析与应用.北京：清华大学出版社，2004附 录求解程序：使用Matlab按照并行算法的思想求解Toeplitz矩阵线性方程组的近似解并分析了它

25、们的误差，充分证明该方法的现实意义和广阔的发展前景。我们只附求三对角非对称的程序如下：n=input（请输入矩阵的维数）；a1=input（请输入主对角线元素）；b1=input（请输入第一行第二列的元素）；c1=input（请输入第二行第一列的元素）；b=zeros(n，1)；b9=1；for i=1:n v=i b(i，1)=b9；endb；A=zeros(n)；for i=1:nfor j=1:n if i=j A(i，j)=a1； elseif j=i+1 A(i，j)=b1； else A(i，j)=0； end end endA；A1=zeros(n)；A1=A；if c1>

26、=b1A1=A/b1；b=b/b1；else A1=A/c1； b=b/c1；endA1；b；a2=A1(2，2)；c2=A1(2，1)；S=solve(f-d=a2，-f*d=c2，f，d)f=S.f；f=double(S.f)；d=S.d；d=double(S.d)；/求，其实由关系我们可以用求解/if abs(d)<1 /使用取绝对值大于1的/ f /由关系，求出/ delse d1=input(您输入元素不符合条件，请重新输入)；endA2=zeros(n)；A2=A1；A2(1，1)=f；A2；m=ceil(n/2) pl=zeros(m)；for i=1:mfor j=1:m

27、 if i=j p1(i，j)=a2；elseif i=j+1 p1(i，j)=c2；elseif j=i+1 p1(i，j)=1；else p1(i，j)=1；p1(1，1)=A2(1，1)； p1(m，m)=A2(n，n)； endend p1k=n-m for i=1:kfor j=1:k if i=jp2(i，j)=a2；elseif i=j+1p2(i，j)=c2；elseif j=i+1 p2(i，j)=1；else p2(i，j)=0；p2(1，1)=A2(1，1)； p2(k，k)=A2(n，n)； end p2s=zeros(n)；for i=1:m forj=1:ms(j，i)=p1(j，i)； endendfor i=1:k； for j=1:k s(m+j，m+i)=p2(j，i)； endendb1=zeros(m，1)；b2=zeros(k，1)；for i=1:m b1(i，1)=b(i，1)；endfor i=1:k b2(i，1)=b(m+i，1)；endb1；b