改进的并行广义共轭残差算法

（，：『“，Ａｐ，）＝（‘・Ａ只）一ｏ：ｓ（Ａｐ，，Ａｐ．）＝ （一一ｌ，Ａ只）一哆一Ｊ（Ａｐ』一Ｉ，ＡＢ）＝
２算法描述
给定线性矩阵方程组
Ａｘ＝ｂ （１）
（（，却，）一ｑ
ＳＡｐ，，Ａｐ．）＝０
（２）
其中，矩阵Ａ∈月”“为非奇异非对称矩阵；ｂ∈Ｒ“，上∈Ｒ”； 针对式（１）给出的ＧＣＲ算法如算法１所示，初始值取
（７）
令ａＪ＝（‘，ｔｊ），ｂｊ＝（ｑ，，ｑｊ），ｃ，＝（ｆ』，ｔｊ），吒＝（ｆ，，ｑ１），
ｉ＝Ｕ＇１，…，』一１，则改进的ＧＣＲ算法如算法２所示。 算法２ ＩＧＣＲ算法
由上述讨论可知，ＧＣＲ算法中的第』步迭代中需要１次 矩阵向量乘、４次向量校正和ｊ＋３次内积计算，ｊ＋１个参数屈， 和１个参数口．校正的浮点计算，该算法的第Ｊ步迭代并行计 算时间为 硭硼＝ｔｍ，＋乞（，＋１）＋‰（２）＋２屯（２）＋２ｔ＂ｆ，＋２，＋（，＋２）ｔ，＝
ｔ。，＝（２ｑ一１）（Ｎ／Ｐ）ｔｒ＋２“，＋ｑｔ。） （５）
ＩＧＣＲ算法总的计算时间为
ｔ，Ｇｃ月：（２ｑ＋３ｎ＋１０）ｆⅣ／Ｐ）ｎｔｒ＋
（２ｎ２十３ｎ）ｔ，十（２ｆ，＋（，ｌ＋５）ｆ。）ｎｌｂＰ＋２ｎ（ｆＴ＋ｑｔ。）
（１ １）
由式ｆ１０）和式（１１）可知，当处理机台数满足下式时
１ｂＰ＞堕竺‘
４ｆ．
（１２）
ＩＧＣＲ算法比ＧＣＲ算法的并行计算时间要少好。
由式（９）求解兰警＝０可得ＧＣＲ算法并行计算时问最小
时的处理机台数
ｋ％群
，｜Ｐ
㈣，
由式（１０）求解ｏｔ、ｔｃ，。ｃ塑＝０可得ＩＧＣＲ算法并行计算时间最
小时的处理机台数 ‰ｃｓ
只。。：—（２ｑ＋３ｎ＋ｌ—Ｏ）Ｎｔ／Ｉｎ ２—百ｉ再矿
２
¨叫 （１４）
在分布式存储并行机中ｔｓ＞＞ｔ。，且２种算法所需迭代次
积，须对这ｋ个局部内积均进行一个全归约（ＡＬＬ—ＲＥＤＵＣＥ） 的通信。其通信时间为２（ｆ。＋ｋｔ。）ｌｂＰ。因此，ｋ个向量内积的 并行计算时间为
（４） ０（≈）＝２（ｋＮ／Ｐ）ｔｒ＋２（ｆ，＋ｋｔ。）ｌｂＰ
（６）
（ｏ，ｑ，＋１）２（％，Ａ乃＋１）＝ｆ‘＂Ａ５。）＋盔岛（０ｌ，Ａｎ）＝（ｒ卅，‰）（３）
（２ｑ＋６ｊ＋１３）（Ⅳ／Ｐ）ｔ，十（Ｊ＋２）ｆ
ｒ＋
（１）ｆｆ算毛＝ｂ—Ａｘｎ，ｔｏ＝Ａｒｏ，０【ｏ＝（ｒｎ，ｔ（Ｉ） 取Ｐｏ＝ｒｏ，ｑｏ＝ｔｏ，ｂｏ＝（ｑｏ，ｑｏ） （２）对ｊ＝０，１，…直到收敛计算；
（３）ａｊ＝０【Ｊ／ｂｊ； （４）ｘｊ＋ｊ
２
２（２ｔ，十（Ｊ＋３）ｔ。）ｌｂＰ＋２（ｔ。＋ｑｔ。）
Ｉｍｐｒｏｖｅｄ Ｐａｒａｌｌｅｌ Ｇｅｎｅｒａｌｉｚｅｄ
（Ｄｅｐａｒｔｍｅｎｔ ｏｆ Ｃｏｍｐｕｔｅｒ Ｓｃｉｅｎｃｅ
Ｃｏｎｊｕｇａｔｅ
Ｐｏｓｔｓ
Ｒｅｓｉｄｕａｌ
Ａｌｇｏｒｉｔｈｍ
ＺＨＡｏ Ｌｉ．ｂｉｎ．ＴＩＡＮ Ｙｏｕ．ｘｉａｎ
ａｎｄ
Ｔｅｃｈｎｏｌｏｇｙ，Ｃｈｏｎｇｑｉｎｇ Ｕｎｉｖｅｒｓｉｔｙ ｏｆ
ｐｅｒｆｏｒｍａｎｃｅ
ｏｆ ｔｈｅ ＩＧＣＲ
ａｌｇｏｒｉｔｈｍ
ｉｓ ｂｅｔｔｅｒ ｔｈａｎ ｔｈａｔ ｏｆ ｔｈｅ ＧＣＲ
ａｌｇｏｒｉｔｈｍ．
［Ｋｅｙ ｗｏｒｄｓ］Ｇｅｎｅｒａｌｉｚｅｄ Ｃｏｎｊｕｇａｔｅ Ｒｅｓｉｄｕａｌ（ＧＣＲ）ａｌｇｏｒｉｔｈｍ；ｐａｒａｌｌｅｌ ｃｏｍｐｕｔａｔｉｏｎ；ｓｙｎｃｈｒｏｎｉｚａｔｉｏｎ
ＸＪ＋ＯｒＪＰＪ；
（８）
ＩＧＣＲ算法中的第．『步迭代中则需要１次矩阵向量乘、
（５）０ｌ ２‘一０【Ｊｑｊ；
（６）ｔＨ ２Ａｔ．１；
４次向量校正和ｊ＋３次向量内积计算，ｊ＋１个参数∥，和１个参 数口，校正的浮点计算，６。校正的浮点计算，该算法一步迭 代并行计算时间为 ｆ南ｃ膏＝０，＋０（』＋３）＋２ｆ＂（２）＋２ｔ坩（』十２）＋（４Ｉ，＋５）ｔｆ＝
ｐｒｏｐｏｓｅｄ ｆｏｌ ｓｏｌｖｉｎｇ
ｔｈｅ ｔｉｍｅｓ
ｎｏｎ—ｓｙｍｍｅｔｒｉｃ ｓｐａｒｓｅ
ｌｉｎｅａｒ
ｃｏｎｖｅｒｇｅｎｃｅ ｏｆ ｃｏｍｐｕｔｅｓ
ＩＧＣＲ ａｌｇｏｒｉｔｈｍ
ｔｈｅ
ＧＣＲ ａｌｇｏｒｉｔｈｍ．ｂｕｔ
ｍｅｍｏｒｙ
ｏｆ
ｔｈｅ ｓｙｎｃｈｒｏｎｉｚａｔｉｏｎ ｏｖｅｒｈｅａｄ
ｏｎ
ａｒｅ ｒ℃ｄｕｃｅｄ
ｂｙ
ａ
ｆａｃｔｏｒ ｏｆ
ｉｔ
ｕｓｉｎｇ
ＩＧＣＲ ａｌｇｏｒｉｔｈｍ
ｄｉｓｔｒｉｂｕｔｅｄ
ｃｌｕｓｔｅｒ ｓｙｓｔｅｍｓ
ｂａｓｅｄ
ＭＰＩ ｅｎｖｉｒｏｎｍｅｎｔ．Ｔｈｅ ｎｕｍｅｒｉｃａｌ
ｒｅｓｕｌｔ ａｎｄ ｔｈｅｏｒｅｔｉｃａｌ ａｎａｌｙｓｉｓ ｐｒｏｖｅ ｔｈａｔ ｔｈｅ
４数值示例
本文数值计算的测试环境为包括８个节点的工作站机 群，主节点的基本配置为Ｉｎｔｅｒ 内存。各节点通过１
０００ Ｐ４ ２．０ ＧＨｚ，５１２ ＭＢ
ｏｖｅｒｈｅａｄ
１概述
天气预报、结构分析、热传导、电力网格分析、地球物 理勘探等科学与工程计算中的许多问题都可归结为求解线性 方程组问题。Ｋｒｙｌｏｖ子空间方法常用来求解稀疏线性方程组， 是许多应用软件的核心算法之一。研究人员对其并行实现已 经作了很多工作。为减少同步开销，文献［１】得到改进的ＣＧ 算法；文献【２］和文献【３］分别讨论如何减少ＣＧ算法在并行机 上的同步开销；文献【４】给出ＣＲ算法的改进形式ＩＣＲ算法。 广义共轭残差算法（ＧＣＲ）１５１是一种求解大型非对称稀疏线性 方程组的Ｋｒｙｌｏｖ子空间方法，它需要向量内积、矩阵向量乘 积和向量校正的计算。在分布式存储的计算机上，矩阵和向 量分布在各台处理机上，向量内积和矩阵向量乘积计算需要 各处理机间的通信，因此，必须进行全局通信的内积计算成 为ＧＣＲ算法并行计算的瓶颈。 本文基于文献【ｌ一４］改变算法计算次序的思想，利用ＧＣＲ 算法的固有性质，针对非对称稀疏线性方程组的求解提出一 种改进的广义共轭残差算法（ＩＧＣＲ）。该算法具有与ＧＣＲ算法 相同的收敛速度，而同步开销次数减少为ＧＣＲ算法的一半。 同时本文在基于ＭＰＩ的分布式存储并行机群上给出数值 示例。
（１０）ＰＨ＝‘州十∑ｋＢ…Ｐ； （１１）ｑＨ＝ｔＨＩ＋∑岛ｐ㈣ｑ；
（１２）结束计算
由于ＩＧＣＲ算法和ＧＣＲ算法具有相同的收敛性，若２种 算法均经过＾次迭代收敛，则由式（７）和式（８）得到ＧＣＲ算法 总的计算时间为
ｔｔ２ｃｅ＝（２ｑ＋３忍＋ＩＯ）（Ｎ／Ｐ）ｎｔ，＋
由算法２可见，内积计算只需要一次全局同步开销。
Ｃｏｎｊｕｇａｔｅ Ｒｅｓｉｄｕａｌ（ＧＣＲ）ａｌｇｏｒｉｔｈｍ
ＧＣＲ
ａｌｇｏｒｉｔｈｍ ｉｓ
ｉｎｔｅｒｄｅｐｅｎｄｅｎｃｅ
ｆｏｒ
ｃｏｍｐｕｔａｔｉｏｎ ｉｎ ｐａｐｅｒ．Ｔｈｅ
ｔｗｏ ｗｈｅｎ
ｔｈｅ ＧＣＲ ａｌｇｏｒｉｔｈｍ，ａｎ
ｉｍｐｒｏｖｅｄ ｐａｒａｌｌｅｌ ｉｓ ｔｈｅ ｓａｍｅ
ｏｎ ａｓ
数较少时有％伸ｚ ２屹曲，可见ＩＧＣＲ算法具有更好的扩
展性。 当ｔ。＞＞ｔ。＞ｔ ｒ，且２种算法所需迭代次数较少时，与 ＧＣＲ算法相比，ＩＧＣＲ算法的性能提高比率为
刀：—ｔＧＣｅ－－—ｔ／ｃ，ｃＲ。
‘ｔＧ曲
！！！！！！￡
（２ｑ＋３ｎ＋１０）Ｎｔ，＋４ｔ。ＰｌｂＰ＋２ｔ。Ｐ
（１５）
对于ｋ个向量内积并行计算时问，每白处理机并行求解 局郎向量内积，计算时间为２（ｋＮ／Ｐ）ｔ，。为计算ｔ个整体内
ｔｈｅ
ａｎｄ Ｔｅｌｅｃｏｍｍｕｎｉｃａｔｉｏｎｓ，Ｃｈｏｎｇｑｉｎｇ ４０００６５）
ａｎｄ ｅｌｉｍｉｎａｔｉｎｇ ｄａｔａ ｌａｒｇｅ
［Ａｂｓｔｒａｃｔ］Ｂｙ ｒｅｌｙｉｎｇ
ｉｎｎｅｒ ｐｒｏｄｕｃｔ ｓｙｓｔｅｍｓ ｉｎ ｔｈｉｓ
ｏｎ
ａｆｔ
ｉｎｔｒｉｎｓｉｃ ｐｒｏｐｅｒｔｙ ｏｆ
Ｇｅｎｅｒａｌｉｚｅｄ
（ｎ２＋３ｎ）ｔ，／２＋（４ｆ，＋（＿ｒｚ＋５）ｆ。）ｎＩｂＰ分析与比较
文献１６ Ｊ对ＧＣＲ算法的收敛性作了详细阐述，并得出结 论：当线性方程组系数矩阵Ａ使ｆＡ＋Ａ７）／２为正定矩阵时， ＧＣＲ算法至多经过Ⅳ步迭代便可给出线性方程组的解。其中 Ｊｖ为矩阵Ａ的阶数。可见对于具有特定结构矩阵的线性方程 组，ＧＣＲ算法所需迭代次数较少。 在分布式存储并行机群结构模型中，假没有Ｐ台处理机。 在２个给定的处理器问完成单位字的一次通信时间为启动时 间ｆ。（包括打包、执行选路算法和建立通信界面的时间Ｊ和传 输单位字的时问ｔ。之和。 以按行划分为例，设每个浮点运算时间均为ｆ，，矩阵阶 数ＪＶ能被处理器台数Ｐ整除，假设不存在负载不平衡的问题。 对一般的稀疏矩阵结构并行矩阵向量乘的纯计算时间为 （２ｑ－］）（Ｎ／Ｐ）ｔ，（每行进行了ｇ—１次加法，ｑ次乘法），其中，ｑ 为每行非０元素的个数，其通信时间与系数矩阵结构以及矩 阵划分方式相关。在最坏情况下，按行划分需要的矩阵向量 乘要用到全部向量，而向量分布存储在每个机器上，因此， 须一次全收集（ＡＬＬ—ＧＡＴＨＥＲ）的通信，通信时间为 ２（ｆ，＋ｆ。Ｎ／Ｐ）ＩｂＰ。而对一般的稀疏矩阵结构，只须相邻处理 机之间进行通信即可，其通信时间为２（‘十ｑｔ。），与处理器台 数Ｐ和矩阵规模Ⅳ无关。下文考虑该种情况的矩阵向量乘通 信时间，其并行计算时间为
（ｑ川，ｑ川）＝（ｆ川，ｔｐｌ）＋２∑岛属（ｆ川，吼）十∑岛霹（ｑｉ，吼）＝ Ｏｆ＋Ｉ，ｔｊ＋ｌ，一∑岛名（吼，吼）
消除ＧＣＲ算法中２次内积计算中的数据相关性。 从武（４）可见，内积（ｑ，，ｑ，）能通过迭代的方法求出，从而
对于ｋ个向量相加的向量校正，其无须通信，因此其并
行计算时间为
ｔｖｅ（≈）＝２（ｋ—１）（＾，／Ｐ）ｔ，
ｘ０＝０。
基金项目：重庆市科委基金资助项目（ＣＳＴ２００５ＢＢｆｇ）６１） 作者倚介；赵利斌（１９８３一），男，硕上研究生，主研方向：并行数值