大规模稀疏矩阵并行计算

第３０卷第６期　２０１０年６月　计算机应用　Ｊｏｕｒｎａｌ　ｏｆ　Ｃｏｍｐｕｔｅｒ　Ａｐｐｌｉｃａｔｉｏｎｓ　Ｖ０Ｉ＿３０　Ｎｏ．６　Ｊｕｎｅ　２０１０　

文章编号：１００１—９０８１（２０１０）０６—１６７９—０３　

基于并行计算的大规模群体行为建模与仿真方法研究　

孟凡亮　，胡晓峰　，蒋亚群　，禹海全　，徐旭林　

（１．国防大学信息作战与指挥训练教研部，北京１０００９１；　２．军械工程学院导弹工程系，石家庄０５０００３；３．南开大学自动化系，天津３０００７１）　（ｙｕａｎｗｈｅｒｅ＠１６３．ｔｏｍ）　

摘要：为解决当前群体行为模型因规模扩大而导致计算量剧增的问题，采用并行离散事件方法构建了大规模　群体行为模型，利用ＹＨ—ＳＵＰＥ仿真引擎实现了群体行为模型的并行计算。重点介绍了模型中仿真对象和仿真对象信　息交互的设计方法，并对该模型在不同数量的节点和仿真实体的环境下进行了测试。实验结果表明，将并行计算引　入群体行为建模之中，可以显著提高仿真个体的数量，更加有效地支持了群体模型的实时运行。　

关键词：并行计算；大规模群体行为；ＹＨ—ＳＵＰＥ仿真引擎；并行离散事件仿真　中图分类号：ＴＰ３９１．９　文献标志码：Ａ　

Ｒｅｓｅａｒｃｈ　ｏｆ　ｍｏｄｅｌｉｎｇ　ａｎｄ　ｓｉｍｕｌａｔｉｏｎ　ｍｅｔｈｏｄ　ｆ０ｒ　

ｌａｒｇｅ　ｓｃａｌｅ　ｃｒｏｗｄｓ　ｂｅｈａｖｉｏｒ　ｂａｓｅｄ　ｏｎ　ｐａｒａｌｌｅｌ　ｃｏｍｐｕｔｉｎｇ　

ＭＥＮＧ　Ｆａｎ．１ｉａｎｇ　，ＨＵ　Ｘｉａｏ．ｆｅｎｇ　，ＪＩＡＮＧ　Ｙａ－ｑｕｎ　，ＹＵ　Ｈａｉ－ｑｕａｎ　一，ＸＵ．Ｘｕ—ｌｉｎ　’　（１．Ｄｅｐａｒｔｍｅｎｔ　ｏｆ　Ｉｎｆｏｒｍａｔｉｏｎ　Ｏｐｅｒａｔｉｏｎ　ａｎｄ　Ｃｏｍｍａｎｄ　Ｔｒａｉｎｉｎｇ，Ｎａｔｉｏｎａｌ　Ｄ咖　ｅ　Ｕｎｉｖｅｒｓｉｔｙ，Ｂｅｌｌｉｎｇ　１０００９１，Ｃｈｉｎａ　２．Ｄｅｐａｒｔｍｅｎｔ　ｏｆ　Ｍｉｓｓｉｌｅ　Ｅｎｇｉｎｅｅｒｉｎｇ，Ｏｒｄｎａｎｃｅ　Ｅｎｇｉｎｅｅｒｉｎｇ　Ｃｏｌｌｅｇｅ，Ｓｈｉｆｉａｚｈｕａｎｇ　Ｈｅｂｅｉ　０５０００３，Ｃｈｉｎａ；　３．Ｄｅｐａｒｔｍｅｎｔ　ｏｆＡｕｔｏｍａｔｉｏｎ，Ｎａｎｋａｉ　Ｔｉａｎｊｉｎ　３０００７１，Ｃｈｉｎａ）　

计算机时代２０１７年第９期　・　３３　・　

ＤＯＩ：１０．１６６４４／ｊ．ｃｎｋｉ．ｃｎ３３・１０９４／ｔｐ．２０１７．０９．０１１　

并行计算在矩阵中的应用　

王涛，赵映诚，刘鑫源　

（中国石油大学（华东）理学院，山东青岛２６６５８０）　

摘要：在解决许多实际问题时，经常需要计算一些高阶矩阵。然而传统的串行计算方法往往效率比较低。因此，需将　

串行程序并行化来提高计算效率。文章分别研究了Ｗｉｎｄｏｗｓ　ＡＰＩ、ＯｐｅｎＭＰ、ＭＰＩ、ＰＰＬ这四种并行计算方法在矩阵乘法　并行化中的应用。通过测试不同规模的矩阵，根据加速比衡量并行化的加速效果，对这四种并行化方法的加速效果进行　

了对比。结果表明，这四种方法都可以提高计算效率，其中ＭＰＩ的加速效果最好。　关键词：计算效率；串行计算；并行化；加速比　

中图分类号：ＴＰ３３８．６　文献标志码：Ａ　文章编号：１００６—８２２８（２０１７）０９—３３—０４　

Ａｐｐｌｉｃａｔｉｏｎ　ｏｆ　ｐａｒａｌｌｅｌ　ｃｏｍｐｕｔｉｎｇ　ｉｎ　ｍａｔｒｉｃｅｓ　

Ｗａｎｇ　Ｔａｏ，Ｚｈａｏ　Ｙｉｎｇｃｈｅｎｇ，Ｌｉｕ　Ｘｉｎｙｕａｎ　

（Ｃｈｉｎａ　Ｕｎｉｖｅｒｓｉｔｙ　ｏｆ　Ｐｅｔｒｏｌｅｕｍ（Ｅａｓｔ　Ｃｈｉｎａ）Ｓｃｉｅｎｃｅ　Ｃｏｌｌｅｇｅ，Ｑｉｎｇｄａｏ，Ｓｈａｎｄｏｎｇ　２６６５８０，Ｃｈｉｎａ）　

Ａｂｓｔｒａｃｔ：Ｉｎ　ｓｏｌｖｉｎｇ　ｍａｎｙ　ｐｒａｃｔｉｃａｌ　ｐｒｏｂｌｅｍｓ，ｓｏｍｅ　ｈｉｇｈ—ｏｒｄｅｒ　ｍａｔｒｉｃｅｓ　ｏｆｔｅｎ　ｎｅｅｄ　ｔｏ　ｂｅ　ｃａｌｃｕｌａｔｅｄ．Ｈｏｗｅｖｅｒ，ｔｈｅ　ｔｒａｄｉｔｉｏｎａｌ　ｓｅｒｉａｌ　ｃｏｍｐｕｔｉｎｇ　ｍｅｔｈｏｄｓ　ａｒｅ　ｏｆｔｅｎ　ｉｎｅｆｆｉｃｉｅｍ．Ｔｈｅｒｅｆｏｒｅ，ｓｅｒｉａｌ　ｐｒｏｇｒａｍｓ　ｎｅｅｄ　ｔｏ　ｂｅ　ｐａｒａｌｌｅｌｉｚｅｄ　ｔｏ　ｉｍｐｒｏｖｅ　ｃｏｍｐｕｔａｔｉｏｎａｌ　ｅｆｆｉｃｉｅｎｃｙ．Ｉｎ　ｔｈｉｓ　ｐａｐｅｒ，ｆｏｕｒ　ｐａｒａｌｌｅｌ　ｃｏｍｐｕｔｉｎｇ　ｍｅｔｈｏｄｓ，Ｗｉｎｄｏｗｓ，ＡＰＩ，ＯｐｅｎＭＰ，ＭＰＩ　ａｎｄ　ＰＰＬ，ａｒｅ　ｓｔｕｄｉｅｄ　ｏｎ　ｔｈｅｉｒ　ａｐｐｌｉｃａｔｉｏｎ　ｉｎ　ｍａｔｒｉｘ　ｍｕｌｔｉｐｌｉｃａｔｉｏｎ．Ｂｙ　ｔｅｓｔｉｎｇ　ｔｈｅ　ｍａｔｒｉｃｅｓ　ｏｆ　ｄｉｆｆｅｒｅｎｔ　ｓｉｚｅ，ａｎｄ　ｍｅａｓｕｒｉｎｇ　ｔｈｅ　ａｃｃｅｌｅｒａｔｉｏｎ　ｅｆｆｅｃｔ　ｏｆ　ｐａｒａｌｌｅｌｉｚａｔｉｏｎ　ａｃｃｏｒｄｉｎｇ　ｔｏ　ｔｈｅ　ａｃｃｅｌｅｒａｔｉｏｎ　ｒａｔｉｏ．ｔｈｅ　ａｃｃｅｌｅｒａｔｉｏｎ　ｅ岱。ｃｔ　ｉｓ　ｃｏｍｐａｒｅｄ．Ｔｈｅ　ｒｅｓｕｌｔｓ　ｓｈｏｗ　ｔｈａｔ　ａｌｌ　ｔｈｅ　ｆｏｕｒ　ｍｅｔｈｏｄｓ　ｃａｎ　ｉｍｐｒｏｖｅ　ｔｈｅ　ｃｏｍｐｕｔａｔｉｏｎａｌ　ｅｆｆｉｃｉｅｎｃｙ，ａｎｄ　ｔｈｅ　ａｃｃｅｌｅｒａｔｉｏｎ　ｅｆｒｅｃｔ　ｏｆ　ＭＰＩ　ｉｓ　ｔｈｅ　ｂｅｓｔ．　Ｋｅｙ　ｗｏｒｄｓ：ｃｏｍｐｕｔａｔｉｏｎａｌ　ｅｆｆｉｃｉｅｎｃｙ；ｓｅｒｉａｌ　ｃｏｍｐｕｔｉｎｇ；ｐａｒａｌｌｅｌｉｚａｔｉｏｎ；ａｃｃｅｌｅｒａｔｉｏｎ　ｒａｔｉｏ　

数值计算功能

向量及其运算

1、向量生成

（1）、直接输入

向量元素用“[ ]”括起来，用空格或逗号生成行向量，用分号生成列向量

a1=[11 14 17 18]

a2=[11,14,17,18]

a2=[11;14;17;18] %列向量

用“’”可以进行向量转置

a1=[11 14 17 18]

a4=a1' %a1行向量，a4列向量

也可以用组合方法：

A=[1 2 3];

B=[7 8 9];

C=[A 4 ones(1,2) B]

（2）、等差元素向量生成

冒号生成法：Vec=Vec0:n:Vecn，其中Vec表示生成的向量，Vec0表示第一个元素，n表示步长，Vecn表示最后一个元素

使用linespace函数：Vec=linespace(Vec0,n,Vecn)，其中Vec表示生成的向量，Vec0表示第一个元素，n表示生成向量元素个数（默认n=100），Vecn表示最后一个元素

vec1=10:5:50

vec2=50:-5:10

vec3=linspace(10,50,6)

2、向量的基本运算

（1）、向量与数的四则运算

向量中每个元素与数的加减乘除运算（除法运算时，向量只能作为被除数，数只能作为除数）

vec1=linspace(10,50,6)

vec1+100

vec2=logspace(0,10,6) %对数等分向量

vec2/100

（2）、向量与向量之间的加减运算

向量中的每个元素与另一个向量中相对应的元素的加减运算

vec1=linspace(10,50,6)

vec2=logspace(0,2,6)

vec3=vec1+vec2

（3）、点积、叉积和混合机

点积：dot函数，注意向量维数的一致性

x1=[11 22 33 44]

x2=[1 2 3 4] a=dot(x1,x2)

sum(x1.*x2) %还可以采用sum函数计算向量的点积

叉积：cross函数，注意向量维数的一致性（由几何意义可知，向量维数只能为3）

[矩阵计算]Lanczos⽅法：求稀疏矩阵特征值

更新： 29 JUL 2016

由知，求矩阵A的特征值，⼤多需要先将其三对⾓化（详细⽅法见徐树⽅先⽣的教材。此处外链），即

T=QTAQ

即找到正交矩阵Q使得T成为三对⾓矩阵。然⽽若A为⼤型稀疏矩阵，常⽤的⽅法如Householder和Givens变换都⽆法充分利⽤A的稀疏性，因此考虑

直接计算T和Q的矩阵元以利⽤A的稀疏性加速运算。

⼀、Lanczos⽅法基本原理

将以上分解式中的Q写成

Q=[q

1,q

2,⋯,q

n]

其中qi为Q的列向量。T写成

T=α

1β

10

β

1α

2⋱

⋱⋱⋱

⋱α

n−1β

n−1

0β

n−1α

n

⽐较

AQ=QT

两边矩阵的每⼀列，得到

Aq

i=β

i−1q

i−1+a

iq

i+β

iq

i+1,i=1,2,⋯,n

由于β

0,β

n,q

0,q

n+1未定义，补充定义β

0q

0=β

nq

n+1=0，这样上式两边乘qT

i得到

α

i=qT

iAq

i,β

i=qT

i+1Aq

i=||Aq

i−β

i−1q

i−1−α

iq

i||

2

可以看出只要给定任意的q

1∈Rn且||q

1||

2=1，就能够利⽤递推关系得到全部q

i,α

i,β

i，迭代格式为

α

1=qT

1Aq

1, //初始值

r

i=Aq

i−α

iq

i−β

i−1q

i−1,

β

i=||r

i||

2, //由上⼀轮α，q产⽣新的值

q

i+1=r

i/β

1 (β

i≠0)

α

i+1=qT

i+1Aq

i+1, i=1,2,⋯,n−1

此即Lanczos迭代。其中qi称为Lanczos向量。在第j步产⽣的矩阵Tj称为j阶Lanczos矩阵，其特征值可能是A的部分特征值的很好的近似。详细内容参

考Krylov⼦空间。

注意其中若某⼀步βi=0，则此时得到的T

j的特征值将都是A的特征值。

⼆、具体算法

Lanczos算法求⼤型稀疏矩阵A特征值（三对⾓矩阵）

1. 输⼊A∈SRn×n,q

1∈Rn (||q

1||

2=1)

2. u

1:=Aq

1, j:=1

3. a

j:=qT

ju

j, r

j:=u

j−a

jq

j, β