（计算数学专业论文）上双对角阵的moorepenrose逆的并行计算.pdf

中 文摘 要 摘要 在这篇文章中，我给出了 三种并行计算上双对角阵矩阵的mo o r e - p e n r o s e 逆的并行算 法。它们是零空间并行方法， g r e v i l le 并行算法和分而治之并行算法。在每一章的结尾， 都给出一个数值例子和一个有关并行效率的定理。同时我还给出了在mp i 环境下上述算 法的源代码。 关键词:上双对角矩阵，并行计算， mo o r e - p e n ro s e 逆，零空间， g re v i l l e 并行算法，分而 治之，mp i 环境 f m a o z h o n g s o h u . c o m第i 页， 共3 2 页上双对角阵的 m- p 逆的并行计算 ab s t r a c t t h r e e p a r a l l e l a lg o r it h m s a r e p r e s e n t e d t o c o m p u t e m o o re - p e n r o s e i n v e r s e o f a b i d i a g o n a l m a t r i x .t h e y a r e n u ll s p a c e p a r a l l e l a lg o r i t h m， g r e v i l l e p a r a l le l a l g o r it h m a n d d iv i d e a n d c o n q u e r p a r a l l e l a l g o r it h m .i n t h e e n d o f e a c h c h a p t e r ,i a l s o g iv e a n u m e r i c a l e x a m p l e a n d a t h e o r e m w h i c h d e a l s w it h c o m p u t in g p a r a l le l e ft fi c ie n c y . a t t h e s a m e t im e , i a l s o g i v e t h e re s o u r c e c o d e o f t h e s e p a r a l l e l a l g o r it h m s i n mp i c i r c u m s t a n c e . k e y wo r d s : b i d i a g o n a l m a t r i x , p a r a l l e l mo o re- p e n r o s e n u l l s p a c e , g r e v i l l e p a r a l l e l a l g o r i t h m , d i v i d e a n d mp i v 7 0 8 2 8 9 关于学位论文使用授权的说明 本人完全了解上海师范大学有关保留、使用学位论文的规定，即:学校有权保留送 交论文的复印件，允许论文被查阅和借阅:学校可以公布论文的全部或部分内容，可以 采用影印、缩印或其他复制手段保存论文。 ( 保密的论文在解密后应遵守此规定) 作 巍 名 : v a l兰 日期 :加 c . s . z 6 。 签 名 : j 呈 一 日期:z b s s . r . 政 p e. : k 第一章 引言 竿 _ 音己i 省 7f 7 _ _ _ r _ j 1 r 广义逆的概念最早是i .f r e d h o l m 于1 9 0 3 年提出的，他给出了积分算子的广义逆，并 称之为“ 伪逆”.1 9 2 0 年， e . h .m o o r e 首先 提出了 矩阵广义逆，他利用投影矩阵 定义了 矩 阵唯一的广义逆。然而，矩阵广义逆真正得到迅速发展，并在各个领域得到广泛应用是 在上世纪5 0 年代以 后。1 9 5 5 年r .p e n r o s e 证明了 mo o r e 所定义的广义逆是满足四个矩阵方程 的唯一 矩阵 1 ,2 ,4 ， 后人称这个唯 一的广义 逆为 m o o r e - p e n r o s e 逆。 如果 a是一个 、 x 、 复 矩阵， 那么a 的 m o o r e - p e n r o s e 逆满足下列四个方程 a x a=a , x a x=x , ( a x ) =( a x ) , ( x a ) =( x a ) 这个唯一的矩阵设其为x ，记为a t . mo o r e - p e n r o s e 逆在许多领域有着广泛的运用，比如 最小二乘问 题、 优化、 迭代法、 数 值分析等， 可以 参看 1 , 2 , 4 。 而m o o r e - p e n r o s e 逆 的计算则构成了 mo o r e - p e n r o s e 逆研究的重要内容。现在比较重要的算法有g r e v i l l e 算法， 满秩分解算法， s v d算法，c l in e 算法，等等。这些算法既有一定的优点，又有一定的 局限 性， 适用范围也不尽相同 可参阅 【 9 , 2 3 。随着现代科学技术的发展， 数值代数计 算的问题也朝着大型化方向发展。国内外己 有许多文献探讨并行计算m o o r e - p e n r o s e 逆的 理论和方法。 1 9 8 8 年， o .y d e v e l a n d e .v k r i s h a n m u r t h y 就在i n f o r m a t i o n p r o c e s s i n g l e tt e r s 2 6 ( 1 9 8 7 - 8 8 ) 发表了并行计算广义逆的论文，具体可参阅 5 。他们的思想是 g r e v i l l e 算 法的并行化。另一些论文处理这个问题的基本思想是首先计算相关矩阵的奇异值分 解，然后对非零奇异值取倒数，最终求一出该矩阵的 mo o r e - p e n r o s e 逆。具体可参看文 献【 3 , 1 0 . 1 1 , 1 2 , 1 3 , 1 4 . 1 5 . 国内早在1 9 8 7 年和1 9 8 8 年，王国荣和陆森泉 就发 表了 并 行 计 算 广 义 逆 a t 和 a 公 n的论 文 【 2 0 . 2 1 a 1 9 9 0 年， 王国 荣 又 发 表了 并 行 计算上 双 对角 阵 a的 广义 逆 a t 的 论文 2 2 = 1 9 9 2 年王国 荣在 i n f o r m a t i o n p r o c e s s in g l e tt e r s 4 1 ( 1 9 9 2 ) 又发 表了 并行计算广义 逆 a t 的 一个改善的 算法 7 1 . 1 9 % 年， 王国 荣 和魏益 民发表了并行计算( m一n ) s v d的论文 8 0 2 0 0 0 年， 魏益民等在 a m c 上发表了关于 用 s u c c e s s iv e m a t r ix s q u a r in g 算 法 并 行 计 算 矩 阵 的a 寿 ， 的 论 文 1 6 a 2 0 0 2 年， 魏 益 民 和 王 国 荣 在 a m c 上 发 表 了 关 于 用 p c r ( 并 行 c r a m e r m 则 ) 计 算 矩 阵 的考丢 的 论 文 1 7 。 本 文 将进一步探讨上双对角阵mo o r e - p e n r o s e 逆的并行计算。 定理1 . 1 设 a e c ， 且 。)n 。 那么 存在正交阵 u b ( 二x 二 ) 和v b ( n x n ) 使 得 (l) u , - a l l , f , 0 d 2 f 2 0 0 d 1 0 dlo00 r.胜resessel;lesl -一 f m a o z h o n g s o h u .c o m第1 页， 共3 2 页上双对角阵的 m- p 逆的并行计算 证明 可参阅 6 1 。 有必 要说明 的是 u b 和v b 都可由 h o u s h o l d e r 矩阵 乘积给出 ， 在 4 中 有 具体的算法。 定理1 .2设a e c ， 且 二)。 。 根据定理 ( 1 . 1 ) 存在正交阵 u b ( m x m ) 和v b ( 二 x n ) 使 得u t a 冷 二t其中 t 是上双对角阵。 那么 a t = v b t t 蜡( 1 . 2 ) 证明 可参阅 1 9 ) 0 由定理1 . 1 和定理1 .2 可知，我们可以对一般的矩阵a e c ，通过正交变换变成 一个上双对角阵t ，这样我们就可以先计算t 的mo o r e - p e n r o s e 逆，然后利用定理1 . 2 计 算 a 的 m o o r e - p e n r o s e 逆。在接下来的第二章将讨论通过计算零空间的方法并行计算奇 异上双对角阵的m o o r e - p e n r o s e 逆，在第三章将用传统的g r e v i l l e 算法并行计算奇异上双 对角阵的 mo o r e - p e n r o s e 逆，以与第一种方法进行比较。当然要对传统的g re v i l l e 算法进 行改造，使其有相当强的并行性。在第四章将讨论用分而治之的方法计算上双对角阵 的 m o o r e - p e n r o s e 逆，以 期得到更 好的并 行性。 本论文中所有算法都有相应的程序并且在s g i 3 0 0 并行机上实验过。程序是在mp i 环境下 编写的，捆绑的语言是c 语言。 具体可参阅 1 8 1 . f m a o z h o n g s o h u . c o m第2 页 ， 共3 2 页上双对角阵的 m - p 逆的并 行计算 第二章 并行计算上双对角阵的 mo o r e - p e n r o s e 逆的零空间方法 第二章并行计算上双对角阵 的mo o r e - p e n r o s e 逆的零空间方法 设 a 为奇异上双对角阵， 上对角元 c , 54 0 ( i 二1 , 2 ， 二 ， n 一1 ) 若 某个 c ; =0 则 把a分 块成下列形式 (zl) 这里a l , a : 均为上双对角阵，并且他们的上对角元均非零，由 1 可知 ( 2 . 2 ) o越 从0 -一 a 故原问题可分块求解。因此不失一般性， 本文假设所讨论的上双对角阵的上对角元均非 零。 下列引 理是【 2 2 中已 证明 的结果。 引理2 . 1设a e r l , r a n k a = r , u, v e r 0 n - r 且r a n k u = r a n k v=二 一: , u 和v 的 列分 别构成 零空间 n ( a t ) 和n ( a ) 的 基 底， 设 (z+3)网 c=a+uv t 则c 非异，且 a t 二c - 1 对于 ( 2 . 1 )中给出的上双对角阵a, 引 理2 .2设上观 对角 奇异阵 一 v ( v t v ) 一 ( u t u 、 一 l tt 容易 求出 n ( a t ) 和n ( 川的 基 底。 r r x n 尹0 “=1 , ， 。 一1 )( 2 . 1 a ) 任 ，.1.esl卫lee.，.j 气-1殊 - 1 b t二 .n r一，.万.月l 一- a f m a o z h o n g s o h u . c o m 第3 页， 共3 2 页上双对角阵的 m- p 逆的并行计算 (习 一 1 ) b l l c l 2 b 1 如 i ( c l c 2 ) ( 一 1 ) n - l b lb_1/(c, b n / ( c l ， . b n /(2 . 侃 _1 n - 1 伪 -2 b 3 以 -l c n -l ( 2 . 5 ) ，.月.1.t!eseej 、t且矛、!廿苦 ( 一 1 ) b , 二 ，b . / 一 : ) 分别为 n ( a ) 和 n ( a t ) 的基底。 证明由 设可知 r a n k a = r a n k a t =。 一1 ， 因此 a x =0 和a t y =0 有非零解x 和 y ，分别构成 n ( a )和 n ( a t ) 的基底。事实上， ( 2 .6 ) 取x 和y 形如( 2 . 1 ) 式， 则 ( 2 .6 ) 式成立。 今把( 2 a a ) 的 a改写成下列分块形式 盯sasb (归但 飞.j 22刀 aa 肉0 rll矛.we.j 一一 a 其中 ，.11 1上门曰;。0 r.r.es.1在才es 11 队口 一一 1 11 a b . - 1 c- l res.1.几l 一- 2 1几 a a 22 一 。 l 0 ，二 ，0 ,1 f m a o z h o n g s o h u .c o m第4 页， 共3 2 页 ( 2 .8 c ) 上双对角阵的m- p 逆的并行计算 第二章 并行计算上双对角阵的mo o r e - p e n r o s e 逆的零空间方法 则 a 括是 。 一 1 阶下 三角阵 ( 2 . 9 ) 。1- n门 00-匀 01一q a i z1 一 c s c a t勺l匀 论一自处勺 ( - 1 ) - z b 3 . .6 。 一 :( - 1 ) - 3 b 3 . .b 。 一 :( 一 w- 4 b , . 一 : b -1 c i c 空 “沈”_1已 2 c - . c .一1c 3 e s. . . c ”一i c -2 c _, 因此( 2 .5 ) 中的 x 和y 可写成下列形式 1 - 6 1 -1 1 z 一 l e r ( 2 . 1 0 x ) y二 - b . a 1 z t e n - 1 1 ( 2 . 1 0 b ) ，1.，les，.eseses 00一1 尸，.三.，.t，lwe刁es 这里 a n d e . , - 1 = 是 。 一1维向 量。 引理2 .3设a 是 ( 2 . 1 a )中给出的上双对角阵， ( 2 . 1 0 a )和 ( 2 . 1 0 b )中给出的x 和y 分别 是n ( a ) 和 n ( a t ) 的 基底， 则 c二 a + y x t件1 1 ) 非异，_目 这里 c - ， 二 d + a /j ( 1 + x d y ) _ y t 一 a 2 2 t d一 o d y y t( 2 . 1 2 ) vj ，.j ud 0一叹 ia d 二 r n x n ， a = ( x t x ) - 1 , q =( y , y ) - :( 2 - 1 4 ) 证明 可参阅 2 2 e 定理2 . 1设a 是、 ( 2 . 1 a )中给出的上双对角阵，则 a t = ( i 一 a x x t ) d ( i 一 q y y , ) ( 2 . 1 5 ) 这里 x 和 y 分 别 是 n ( a t ) 和n ( a ) 的 基底， 由 ( 2 .5 ) 给出 ，d , a , ,3由 ( 2 . 1 4 ) 给出 。 证明 由引理2 . 1 - 2 . 3 知 a t =d 十 。 侧 1 + x t d 刃 x y ， 一 二x t d一 o d y y t 一 a o x 沪 f m a o z h o n g s o h u .c o m第5 页 ， 共 3 2 页上双对角阵的 m - p 逆的 并 行计算 =( i 一 。 x x ) d ( i 一 q y y , ) . 上述定理表明，上双对角阵的m- p 逆具有一个相当简单的结构，且适合于并行计算。 利用( 2 . 1 0 ) ,( 2 . 1 1 ) ,( 2 . 1 4 ) 和( 2 . 1 4 ) ， 我们可得到 算法2 . 1设a 为n 阶奇异上双对角阵，使用n 台处理机并行计算 a t 的步骤为 1 ) 并 行 计 算a 扮 . 2 ) 并行计算x和 ， 3 ) 并行计算 a和 尽 4 ) 并 行 计算( i 一 a 二 二 t ) d . 5 ) 并 行计算 ! ( i 一 a x x t ) d ( i 一 q y y , ) =a t . 下面讨论算法2 . 1 的并行复杂性。 定 理 2 .2设 a 是( 2 . 1 a ) 中 给出的 上双对角阵， 使用 n 台 处理 机， 算法 2 . 1 能在 o ( n r l o g n l ) 步内 算出 a t， 速度增长倍数为 民= o(- f l o g n l ) 证明 1 ) 并 行 计 算a 盆. 由( 2 .9 ) 知 ， 我 们 可 以 用 下 列 方 式 逐 次 计 算n 一 1 阶 阵a 衬的 各 条下对角线上的元素。为简单起见，列出n 一1 二4 的情况。 c l c 2 c s c q 1-伪 1-c3 一b 3 1 - b 3 c 2 c 3 1-处 一b 2 1 - b 2 c 1 c 2 一b 4 1 - b 4 c 3 c 4 1一cl 鱼晌 4-c -b.一q伪-冰 业cac3-ba一。 渔、竺clczc3 业c1c2 因此，使用n 一1 个处理器能在c t , 二 - b 2 b 3 b 4 c i c 2 c 3 c 4 。 步内 算出4 1 21 . 2 ) 并 行计算 x 和 y 。因a 扮的 全 部元素均己 算出 ，( 2 . 1 0 )，使用n 一1 个处理器作一步 由= 乘法算出x . 3 ) 并行计算 再作一步乘法可算出 y 。因此，c t z 。和 iq因 n 维向 量x 和 y 的分量x ; 和y i 均己 算出，由( 2 . 1 3 ) 器作一 步乘法和 ( l o g n l步加 法可 算出 艺m 2 ， 再作一步除 法便可算出a : =1 q. 因 此，c t 3 =2 ( f l o g n l + 2 ) . 4 ) 并 行 计 算 ( i = a x x t ) d = g冷 使用n 个处理 同理可算出 a i z1 =( w 1 , w 2 , ， 二 ， - n - 1 ) 因为a 讨为。 一 1 阶下三 角阵， 故岭的 前7 一1 个 分量为 零。 从 而 d 的 翁列 的蒯个分量为零。 利用( i 一 。 ) x x t 的 特点， 与d 的前。 一i 列逐步相乘，使用n 个处理 f m a o z h o n g s o h u .c o m第6 页， 共3 2 页上双对角阵的m - p 逆的并行计算 第二章 并行计算上双对角阵的m o o r e - p e n r o s e 逆的 零空间方法 器，各种计算所需步数如下:a x作一步乘法，作一步乘法和 f l o g ( 二一 a l 作一步减法; 因此，c t q = 、!矛 0哟 了矛!、 步加 法;a x x t作一步乘法; 2 ( f l o g ( n 一 j ) l + 4 ) . 5 ) 并行计 算 g ( i 一 )3 y y t ) =a t . 记 n 阶阵 g 为 ，.lweeseseseswewej tt，t 91乳如 尸卫.11.1.l 一一 g 这里g i t 是 n 维行向 量， 利 用i 一 口 ， 沪的 特点 计算g : 叹 i 一 q y y 勺、 使 用 n 个处理 器， 各 种计算所需步数如下:o y t 作一步乘法:g i t ， 作一步 乘法和 f l o g n l步加法;g ; t y q y t 作一步 乘法;g i t 一 g i t 妇尹作一步减 法。 因 此，c t s =ef l o g - +4 =n ( f l o g n +4 ) i =1 艺c t ; 一 。 + 2 + 2 f lo g e + 4 + 艺( 二 (n 一 j )1 + 4 ) + n ( f lo g n + 4 ) 为估计上式，设 l o g ( 。 一1 ) l =4因为 又flo g ( - 一 a l - 艺f lo g p l p =2 f l o g 2 l +r l o o 3 l + r l o g z 2 1 +一 +r ( l o g 2 - + 1 ) p + 二 + r l , , 2 1 +. +r ( l o g 2 ， 一 ， + 1 ) 1 + 二 +r l o g ( 。 一1 ) 1 + 2 1 * 2 + 2 2 * 3 + + 2 q - 2 * ( q 一 1 ) + ( n 一 1 一 2 q - 1 ) q ( 1 + 2 + 2 2 + + 2 q - 2 ) + ( 2 + 2 2 + + 2 q - 2 ) + + ( 2 q - 3 + 2 q - 2 ) + 2 q - 2 +( 。 一 1 一 2 q - ) q ( 。 一1 ) q 一 2 q +1 . 所以 算法 ( 2 . 1 ) 能 在e c t = o ( r l o g n ) 步内 算出a t . 因 为相应的串 行算法 在单处 理器 上的 运行时间为o ( n 2 ) 故速度增长倍数为 凡=o(- r l o g n l ) f m a o z h o n g s o h u .c o m第7 页，共3 2 页上双对角阵的 m- p 逆的并行计算 数值例子: a 160000 025000 003400 000030 000042 000005 ( 2 . 1 6 ) 用6 个c p u 计算得出其的mo o r e - p e n r o s e 逆为: 0 . 0 3 3 5 5 7 - 0 . 0 2 0 1 3 4 0 . 0 1 2 0 8 1 0 0 0 0 . 1 6 1 0 7 4 0 .0 0 3 3 5 6 - 0 . 0 0 2 0 1 3 0 0 0 -0 . 0 6 4 4 3 0 0 . 1 9 8 6 5 8 0 . 0 0 0 8 0 5 0 0 0 0 . 0 4 8 3 2 2 - 0 . 1 4 8 9 9 3 0 . 2 4 9 3 9 6 0 0 0 0 0 0 0 . 1 3 1 6 1 9 0 . 1 5 1 2 8 6 - 0 . 0 6 0 5 1 4 0 0 0 - 0 . 0 3 6 3 0 9 0 . 0 2 7 2 3 1 0 . 1 8 9 1 0 7 ( 2 . 1 7 ) 用时为 0 . 8 1 5 7 3 7 s e c o n d s . 而用m a t l a b ( 6 .5 ) 计算的结果是: 0 . 0 3 3 6 -0 . 0 2 0 1 0 .0 1 2 1 0 0 0 0 . 1 6 1 1 0 . 0 0 3 4 - 0 . 0 0 2 0 0 0 0 -0 .0 6 4 4 3 0 0 . 1 9 8 6 5 8 0 . 0 0 0 8 0 5 0 0 0 0 . 0 4 8 3 -0 . 1 4 9 0 0 . 2 4 9 4 0 0 0 0 0 0 0 . 1 3 1 6 0 . 1 5 1 3 -0 . 0 6 0 5 0 0 _0 - 0 . 0 3 6 3 0 . 0 2 7 2 0 . 1 8 9 1 ( 2 . 1 8 ) leseseseses.es，.月.es.esllwej 下面是用6 个c p u 计算6 阶上双对角方阵的mo o r e - p e n r o s e 逆的m p i 并行程序: #i n c l u d e m p i . h # i n c l u d e 禅 i n c l u d e i n t m a i n (i n t a r g c , c h a r * a r g v ) 【 i n t m y i d , n u m p r o c s , i , j ; d o u b l e a 5 5 ， t 5 ( 5 ) ; d o u b l e s t 6 6 , g 6 6 1 , m 6 ， f ( 6 : d o u b l e b 6 ， c( 5 , d 5 ; d o u b l e x 6 , y ( 6 , k 6 ， r 6 , w 6 ， h 6 ; d o u b l e p, q; d o u b l e s t a r t t i m e, e n d t i m e; f m a o z h o n g s o h u .c o m第8 页， 共3 2 页上双对角阵的 m - p 逆的并行计算 第二章 并行计算上双对角阵的 mo o r e - p e n r o s e 逆的零空问方法 mpi mpi 从pi i n i t ( s i z e ( m p i r a n k( m p icomm w o r l d , w o r l d , i f ( m y i d = = 0 ) p r i n t f ( p l e a s e i n p u t f o r( i = o ; i 6 ; i + + ) s c a n f ( 11 % l f ， i 5 ; i + + ) s c a n f ( e l f ， p r i n t f ( .i n ) ; v e c t o r b n ) ; v e c t o r c n ) ; m p i $ a r r i e r( m p i c o m m w o r l d) ; s t a r t t i m e = m p i w t i m e( ) ; f o r ( i = 0; i 5 ; i + + ) d ( i = o ; f o r( i = o ; i 6 ; i + + ) ( x i = 0 ; y( i 二 0 : kf i = 0; r i l = o ; w i j = 0 ; h i = 0 ; f o r ( i = 0 ; i 6 ; i + + ) f o r ( j 二 o ; j 6 ， j + + ) s ( i j = 0 ，g i l ( j = 0 ; 。 i 1 j 1 = 0 ; # ( i l j l = o ; f o r( i = o ; i s ; i + + ) f f o r ( j = o ; j s ; j + + ) a i l j l = o ; t i j ) 二 o ; f m a o z h o n g s o h u . c o m第9 页，共3 2 页上a对角阵的m- p 逆的并 行计算 p = 0 q二0 m p i b c a s t( c o m m w o r l d ) ; m p i b a r r i e r( m p i c o m m w o r l d ) ; f o r ( i = 0 ; i s ; i + + ) d i l = o; i f ( m y i d ! = 5 ) d( m y i d l = l / ( c ( m y i d ) ) ; f o r ( j = m y i d + l ; j 5 ， j + + ) d ( j 二 ( 一 1 ) * d ( j 一 1 1 * b( j l / c ( j l ; m p i b a r r i e r( m p i c o m m m p i 一 a l l g a t h e r ( d , 5 , mpi mpi f or b a r r i e r( m p i c o m m w o r l d) ; d o ti 卫l e, w o r l d); ( i = 0 ; i 5 ; i + + ) f o r ( j = o ; j 5 ; j + + ) t i l j = a j i l ; m p i b a r r i e r( m p i c o m m w o r l d ) ; x( 0 1 = 1 ; y 5 = 1 ; f o r ( i = 1 ; i 6 ; i + + ) x i ) = ( - b ( o l ) * t i 一 1 1 0 l ; y ( 1 一 1 1 = ( - b ( 5 1 ) * a i 一 i l ; m p i _ b a r r i e r ( m p i _ c o m m 一 w o r l d ) : p = o ; q = 0 ; f o r( i = o ; i 6 ; i + + ) f m a o z h o n g s o h u .c o m第1 0 页， 共3 2 页上双对角阵的m- p 逆的并行计算 第二章 并行计算上双对角阵的 mo o r e - p e n r o s e 逆的零空间方法 p + = x i * x i ; q + = y i * y i : p = l / p ; q = l / q ; m p i b a r r i e r( m p i c o m m w o r l d) ; f o r( i = 0 ; i 6 ; i + + ) (。 0 i = 0 : 。 i 5 二 0 ; f o r( i = 1 ; i 6 ; i + + ) f o r ( j 二 o ; j 5 ; j + + ) 。 ( i ; j = t i 一 1 3 j ; m p i 一 a r r i e r ( m p i _ c o m m 一 w o r l d ) ; f o r( i = o, i 6; i + + ) k i ! ( - p ) * x m y i d * x i ; k m y i d 二 1 + k m y i d ; f o r( i = 0 ; i 6 ; i + + ) f o r ( j 二 0 ; j 6 ; j + + ) r i + = k j * s j i ; m p i b a r r i e r ( m p i c o i v w o r l d ) ; m p i 一i l g a t h e r ( r , 6 , m p 工 夕o u b l e , m p i b a r r i e r( m p i c o n l m w o r l d ) ; f o r( i = 0 ; i 6 ; i + + ) w ( i 二 ( - q ) * y( m y i d w m y i d = 1 + w m y i d * y ( i ; m p i a l l g a t h e r ( w , 6 , m p i 一o u b l e , m p i b a r r i e r( m p iw o r l d) ; f o r( i = o ; i 6 ; i + + ) f o r ( ? = j ; j 二 ; j + + 1 f m a o z h o n g )s o h u . c o m 第1 1 页洪3 2 页上双对角阵的m - p 逆的并行计算 h i + = g ( m y i d ) j * m j i ; m p i ee a l l g a t h e r ( h , 6 , m p i _ d o u b l e , m p i b a r r i e r( m p i c o m m w o r l d ) ; e n d t i m e = m p i w t i m e( ) ; i f ( m y i d = = 0 ) 丈 f o r( i = 0 ; i 6 ; i + + ) f o r ( j = o ; j 6 ; j + + ) p r i n t f ( k i t， ， f ( i ( j ) ; p r i n t f ( . n . ) ; p r i n t f ( t h a t t o o k s t l f s e c o n d s n ， e n d t i m e 一 s t a r t t i m e ) ; m p i f i n a l i z e( ) ; f m a o z h o n g s o h u . c o m第1 2 页， 共3 2 页上a对角阵的m- p 逆的并行计算 第三章并行g : e v i l l e 递推算法 第三章 并行g r e v i l l e 递推算法 为了与算法2 . 1 作比 较，我们给出了下列并行g r e v i l l e 递推算法。 算法3 . 1设a 为n 阶方阵，使用n 台处理机并行计算 a t的步骤为 1 ) 若。 ; 76 。 ， 并 行 计 算a l t = ( a , t a i ) - l a , t ， 否 则a l t = 0 . 2 ) 对 k =2 , 3 , . . , n. 1 ) 并 行 计 算丸 = a 、 一 户 a k , 2 ) 并 行 计算 c k =a 。 一 a 。 一 ， d k . 3 ) 若c k 7 。 ，并 行 计 算b k t = ( c k t c k ) - 1 c k t , 否 则 并 行 计 算b k t = ( 1 十 d k t d k ) - d k t a * 一 扒 4 ) 并行计算 a k t = d = a k _ i t 一 d k b k 7 b k t 3 ) a n t 二a t 下面讨论算法3 1 的并行复杂性。 定理 3 . 1 设 a 为 ” 阶方阵， 使用 n 台 处理 机， 算法 3 . 1 能 在 o ( n 2 卜 v g 川) 步内 算出川 度增长倍数为 氏=o ( 可l o g n ) 证明1 ) 并 行 计 算a i t ， 要 作 3 步 乘 除 法 和 l o g e 步 加 法， 因 此c t , =f l o g n l + 3 2 ) 并行计算d k , c k , b k t 和a k t ，对每个 k ，并行计算d k 要做k 一1 步乘法和 ( k - 1 ) l o g n l 步 加 法; 并 行 计 算c k 要 做二 步 乘 法 和。 f l o g ( k - 1 ) 十1 步 加法: 并行 计 算 纵 ， 当。 、 = a 时 要 做。 + 3 步 乘 法 和( 二 + 1 ) f l 姆 恤一 川十 1 步 加 法 : 并 行 计 算a 扩 要 做 k 一1 步乘法和k 一1 步加法:因此 c t , =e ( k 一 1 ) ( f l o g n l +1 ) + f l o g ( k 一 1 ) l + 1 ) + 1 + ( n +1 ) l o g ( k 一 1 ) 1 + 二 + k =2 4 + 2 ( k 一 1 ) ， s 因 此算 法 3 . 1 能在 ec t , 二0(n2 l o g 川 ) 步内 算出才 。 因 为 相 应的串 行算 法 在 单处 理 i =1 机 上 的 运 行时 间 为o ( n 3 ) ，故 速 度 增 长 倍 数 为凡= o f 可 l o g e i ) 因此，算法3 . 1 和算法2 1 的速度增长倍数相同，即两个算法的并行性相同。但就计算上双 对角阵之m, 逆而言，算法3 . 1 的运算步数少于算法2 . 1 . f m a o z h o n g s o h u . c o m第1 3 页， 共3 2 页上双对角阵的 m- p 逆的并行计算 160000 025000 003400 000030 000042 000005 r.11，.1. 一一 a 数值例子: 用6 个c p u 计算得出其的m o o r e - p e n r o s e 逆为: 0 .0 3 3 5 5 7 - 0 . 0 2 0 1 3 4 0 . 0 1 2 0 8 1 0 0 0 0 . 1 6 1 0 7 4 0 . 0 0 3 3 5 6 - 0 . 0 0 2 0 1 3 0 0 0 -0 . 0 6 4 4 3 0 0 . 1 9 8 6 5 8 0 . 0 0 0 8 0 5 0 0 0 0 . 0 4 8 3 2 2 - 0 . 1 4 8 9 9 3 0 . 2 4 9 3 9 6 0 0 0 0 0 0 0 . 1 3 1 6 1 9 0 . 1 5 1 2 8 6 - 0 .0 6 0 5 1 4 0 0 0 - 0 . 0 3 6 3 0 9 0 . 0 2 7 2 3 1 0 . 1 8 9 1 0 7 ( 2 . 2 ) ，!一t!les.，j 用时为 0 .0 6 3 5 0 2 s e c o n d s . 而用m a t l a b ( 6 . 5 ) 计算的结果是: 0 . 0 3 3 6 0 . 1 6 1 1 -0 . 0 6 4 4 3 0 0 . 0 4 8 3 0 0 -0 . 0 2 0 1 0 . 0 0 3 4 0 . 1 9 8 6 5 8 -0 - 1 4 9 0 0 0 0 . 0 1 2 1 -0 . 0 0 2 0 0 . 0 0 0 8 0 5 0 . 2 4 9 4 0 0 0 0 0 0 0 . 1 3 1 6 -0 . 0 3 6 3 0 0 0 0 0 . 1 5 1 3 0 . 0 2 7 2 0 0 0 0 -0 . 0 6 0 5 ( ) - 1 8 -1 1 ( 2 . 3 ) 下面是用6 个c p u 计算6 阶方阵的m o o r e - p e n r o s e 逆的m p i 并行程序: # i n c l u d e m p i . h # i n c l u d e # i n c l u d e i n t m a i n (i n t a r g c , c h a r * a r g v ) i n t m y i d , n u m p r o c s , i , j , k ; d o u b l e a 6 ) ( 6 1 , g ( 6 1 ( 6 ) ; d o u b l e e 飞 6 ) 6 ) ; d o u b l e b 6 ) ， c 6 1 ， d 6 ) , m 6 ) ; d o u b l e s , l , q , r , v, w , p a r t s u m ; d o u b l e s t a r t t i m e, e n d t i m e ; m p 工 工 n i t( f m a o z h o n g s o h u .c o m 第1 4 页， 共3 2 页 上双对角阵的 m- p 逆的并行计算 第三章并行g r e v i l l e 递推算法 mpi m pi 三 o m m l : o m m s i z e ( m p i r a n k( m p i comm comm world world i f ( m y i d = = 0 ) p r i n t f ( p l e a s e i n p u t m a t r i x a n ) ; f o r( i = 0; i 6; i + + ) f o r ( j = o ; j 6 ; 千 + ) s c a n f ( e l f ， f o r( i 二 0 ; i 6 ; i + + ) m p i _ b c a s t ( i