（工程力学专业论文）子结构级并行计算的若干问题的研究.pdf

摘要 y - 3 6 2 7 2 2 伪满足大规模科学计算和工程分析的需要，可利用多台微机通过局 、，一1 ， 域网构建的松散耦合多处理机实施并行计算了论文在相关的四个方 面开展了研究，使得有限元并行计算能够更加适合这种特定并行计 算环境。第一，提出了一种以基点为基础的子结构自动分解方法。 艏先由“拓扑距离极大化”的准则确定基点，然后从基点开始采用 抢占方式收集相邻的单元生成长子结构，能灵活控制子结构的数目 和尺度。j 第二，采用静态子结构的凝聚技术，在结构级实行有限元 分析并行化。f 大幅度压缩了节点机之间的数据交换量，提高了有限 7 “7 元并行计算的效率o 第三，利用广义特征值f q 题子空i n 迭代法的算 法特点，结合静态子结构凝聚技术，设计出子空间迭代法的子结构 实现方案，用于结构模态分析可以得到与整体结构完全等价的计算 ，、， 结果。坼且结果精度不受子结构分解形式影响。j 第四，提出用透明 效果表现有限元模型的可视化方案。( 不仅可以显示三维实体中高应 、 力区表面分布，而且可以显示其深度，增强了可视化表现力。卜 关键字：子结构自动分解；并行计算；静态凝聚；子空间迭代；模 态分析； a b s t r a c t i tc a nb ei m p l e m e n t e db yr e l a xc o u p l i n gm u l t i p l ep r o c e s s o r sv i al a nc o n n e c t i o nf o r t h er e q u i r e m e n to fc o s m i c a l l ys c i e n t i f i cc o m p u t a t i o na n de n g i n e e r i n ga n a l y s i s t h e p a p e rw o r k e do v e rf o u r c o r r e l a t i v e a s p e c t s ，w h i c hm a k i n gf i n i t e e l e m e n tp a r a l l e l c o m p u t a t i o nf i t t h es p e c i a lp a r a l l e l c o m p u t a t i o nc i r c u m s t a n c em u c hm o r e f i r s t ，a n e ws u b s t r u c t u r ea u t o m a t i s md e c o m p o s i o nm e t h o db a s e do nr o o tv e r t e xw a sb r o u g h t f b r w a r d ，a n dr o o tv e g e xw a sf i x e do nt h er u l ew h i c hg r a p hd i s t a n c ei sm a x i m u m f i r s t l y , t h e ns u b s t r u c t u r e sw e r em a d eb ya d j a c e n te l e m e n t sc o l l e c t e db yr o o tv e r t e xi n s n a t c h i n gs t y l e ，s ot h em e t h o dc a nn e a t l yc o n t r o ls u b s t r u c t u r ea m o u n ta n dm e a s u r e s e c o n d l y ，s t a t i cs t a t es u b s t r u c t u r ec o n d e n s i o nt e c h n o l o g yw a si n t r o d u c e d ，a n df i n i t e e l e m e n t a n a l y s i sw a s c a r r i e do u t p a r a l l e lo n s t r u c t u r ec l a s s ，s oe x c h a n g ed a t ab e t w e e n m a c h i n en o d e sw a sc o m p r e s s e dg r e a ts c o p ea n dp a r a l l e lc o m p u t a t i o ne f f i c i e n c yw a s i m p r o v e d t h i r d l y ，u s i n gs u b s p a c ei t e r a t i v em e t h o do fg e n e r a l i z e de i g e n v a l u ea n d c o m b i n i n g s t a t i cs t a t es u b s t r u c t u r ec o n d e n s i o n t e c h n o l o g y , t h e s c h e m eb a s e d s u b s p a c e i t e r a t i v ea n ds u b s t r u c t u r em e t h o dw a sw o r k e d o u t b y t h ee x a m p l e so fm o d e a n a l y s i s ，r e s u l ti se q u i v a l e n tt oi n t e g e rf r a m e ，a n dr e s u l tp r e c i s i o nw i l ln o ta f f e c t e db y t h es u b s t r u c t u r ed e c o m p o s i t i o ns t y l e t h i sp a p e rp r e s e n t sav i s u a l i z a t i o nm e t h o do f 3 df i n i t ee l e m e n ta n a l y s i st h e h i g hs t r e s sz o n ei sf i l l e dw i t h t r a n s l u c e n tm a t e r i a la n d h el o ws t r e s sz o n ei sf i l l e dw i t l lt r a n s p a r e n tm a t e r i a l i nt h i sw a y , t h e3 di m a g eo f h i g hs t r e s sz o n e ，b o t ht h es h a p e o ns u r f a c ea n dt h ed e p t ho ft h ez o n e ，c a nb eo b s e r v e d d i r e c t l y k e y w o r d s ：s u b s t r u c t u r ea u t o m a t i o nd e c o m p o s i t i o n ；p a r a l l e lc o m p u t a t i o n ；s t a t i c c o n d e n s a t i o n ；s u h s p a c el t e r a t i v em e t h o d s ；m o d ea n a l y s i s 1 引言 计算机c p u 性能不断提升，i n t e l 和a m d 开发的主频达1 g h z 的芯片即将 面世，不久后具有每秒l o 亿次计算能力的个人计算机将出现在家庭和办公室的 桌面上，十多年以前这是不可想象的。与此同步，计算技术和软件设计方法也 得到了长足发展，性能卓越的c a e 软件日益成熟，不同功能的软件从概念设计、 工程分析、可视化和结构优化实现了无缝连结。由高性能计算机和通用c a e 软 件架构的虚拟试验实验室( v i r t u a le x p e r i m e n t a ll a b o r a t o r y ) 十f l 当程度上可以代替 传统意义上的实验测试。 尽管如此，经验和现实都显示，计算机的性能始终不能满足同时代科学计 算和工程分析的需要。七十年代在运算速度不足百万次7 1 9 计算机上，采用少 于1 0 0 自由度的超参厚壳有限单元模型，分析航空发动机叶片的固有频率已基 本上能够满足设计要求。到了九- 卜年代仞，设计要求分析的单叶片自由度数已 经超过1 0 ，0 0 0 。而且航空发动机结构完整性设计要求分析的对象不仅是根部固 支的叶片，还需要进行盘片、盘片轴，乃至包括转子、静子和机匣的发动机 整体系统的分析。民用汽车车身的模态分析要求也达数十万自由度。高性能计 算始终是工程设计所依赖的工具，长期以来科学和工程计算都是推动计算机发 展的主要动力。 现代设计需要高性能计算机支持，但是所付出的费用也是巨大的。在p c 机迅速发展的今天，与p c 相比高性能计算机的价格显得十分昂贵的，并不是 所有需要大规模计算的企业或研究机构都能承受的。在此背景下，将多台p c 通过网络组合成“松散耦合多处理机”应运而生。由于p c 价格低廉，所以利 用多台p c 构建的多处理机系统开展大规模科学和工程计算成本不高。p c 具有 良好的性能价格比、软件丰富、通用性强，并行计算时作为多处理器系统节点 机使用，平时则可以用于办公自动化、数据管理，或者是教学培训，不会出现 闲置，所以设备的利用非常有效。近年来出现了一批基于消息传递的并行程序 设计的软件平台( 如p v m 、m p i 等) 使得操作由p c 构建的多处理机如同是在u n i x 下操作一台真正的多处理器并行汁算机，并行计算的程序编程和运行也十分便 捷，同时也保证了并行程序的可移植性，能够毫无困难地在同类型并行操作系 统支持下的多处理机之间移植。 当然，由p c 构建的多处理机是一种松散的耦合系统，其特征与真正的并 行计算机是有差异的，其中有些还与算法的设计相关： 夺节点机之间数据传递是通过网络完成的，相对于节点机内部通过总线交 换数据交换速度要低很多。而且局域网上的数据交换是点对点方式进行 的，也就是当两台节点机交换数据时，其它节点机之间的数据交换必须 处于等候状态。所以数据交换是并行计算的瓶颈所在，算法设计应当考 虑减少节点机之间数据交换量。 夺节点机数目是任意的，有别于真正的多处理机中的2 ”格式。如果硬性 满足2 格式要求，就会浪费设备资源。例如原本3 0 台p c 可以作为节 点机使用，为了满足2 ”格式就只有1 6 台p c 能够得到利用，其余1 4 台 只能闲置。所以并行计算应突破2 ”格式的限制，充分发挥设备的效能。 夺用作节点机p c 机群往往是在不同时期购置的，所以无论是c p u 频率、 还是硬盘容量都是有差异的，这也有别于真正的多处理机。因而节点机 负荷平衡不应当作为并行算法的设计要求，而是应当根据节点机的实际 性能合理安排任务，尽可能使得节点机的工作进程保持同步。 针对p c 构建的松散耦合多处理器系统的特性，本文在三方面进行了研究， 以提高有限元素并行分析的有效性。 1 静态子结构并行计算 大规模有限元分析并行化最简单的方法是在平衡方程组求解的层次上实 旆，而且有很多现成的并行算法可以利用。2 f 。但在这一层次实施，对数据交换 的要求很高。如果刚度矩阵的元索为，节点机数量为p ，那么需要传输的数据 量就是w ( p 一1 ) ，而且点对点的传输方式又不能并行进行。 为了减少数据传输量，可以采取静态子结构方法，在结构级实行有限元分 析并行化。策略是： 夺由主机将子结构有限元数据发送到各个从机。 夺各从机并行生成刚度矩阵和载荷向量；实施静态凝聚、将子结构自由度 凝聚到界面上；然后将凝聚后的刚度矩阵和载荷向量发回到主机。 冷由主机进行耦合分析，得到界面位移结果，并再次发送到各个从机。 夺各从机并行计算子结构内部位移和单元应力。 由此可见，主机与从机之间传递的只是有限元模型数据、凝聚后刚度矩阵、 载荷向量和界面位移向量，相对于整体结构的刚度矩阵而言数据量是微不足道 的，容易解决松散耦合多处理系统中存在的数据传输瓶颈。 此外，即使是在串行计算模式下，予结构静态凝聚技术本身也是提高计算 效率的途径。所以静态子结构并行计算方法既能实现有限元分析并行化，又可 以提高有限元分析的效率，尤其适合在松散耦合多处理系统上运行。 2 子结构自动分解 有限元分析的网格应当自动生成，一般可以借助通用c a e 程序。通用c a e 程序得到的往往是最终的生成结果，为了采用子结构并行计算就有必要在此基 础上自动地分解予结构。文献 7 引介绍了若干予结构自动分解算法，但均要求满 2 足如2 m 格式的限制，某些算法还需要进行连结矩阵的特征值计算“”3 。 考虑到p c 构建的多处理机特点，子结构分解应满足： 夺子结构数目与投入运算的节点机一致，提高设备利用率。 夺子结构大小与节点机的性能相匹配，达到物尽其用。 夺子结构界面尽可能小，减少数据传输量。 然而找出“最佳”的分解结果几乎是一个不可解问题，为寻找“最佳”结 果付出的代价可能比有限元分析还要高，这样的“最佳”结果是没有实际意义 的，实际应用中找出一个“较佳”方案也就能够满足了。 分解过程分两阶段进行： 夺在结构边界上确定p ( 节点机个数) 个基点。 令以基点为基础，采用抢占方式收集相邻的单元生长子结构。 确定基点的准则是使得基点之问的最短拓扑距离极大化，收集单元时根据 预先要求的子结构规模设定优先级，保证子结构按比例生长。 3 子空间迭代法的子结构实现 动态子结构方法口，5 1 剐是解决复杂结构模态分析行之有效的方法。这些方法 都假设了子结构已经存在，一般认为子结构就是物理意义上的部件，如飞机是 由机翼、尾翼、垂直翼及若干段机身等自然部件构成的子结构组合体。但如果 子结构是在整个结构基础上采用自动分解算法生成的，那么子结构就不再是自 然的物理部件。经验认为，非自然部件的动态子结构方法会影响计算精度。 在对比各种常用广义特征值问题算法的基础上，发现子空间迭代法可以用予 结构实现。子空间迭代法的计算量集中在两个步骤： 夺逆幂法修正迭代向量 k u k + 1 = m u 夺产生缩聚刚度和质量矩阵 k k + l = u l l k u m = 叫+ 1 m u 逆幂法由d 。修正u 。形式上完全等同：f 在惯性力m o 。作用下求结点位移u 。 的静力学问题，于是前述的“静态子结构并行计算”方法可以将逆幂法并行化。 而缩聚刚度和质量矩阵( 。和m 。与结构的弹性势能和动能有关，可以按子结 构累加计算，所以也能够进行并行计算。 在本质上，子空间迭代法的子结构实现完全等价于子空间迭代法的整体结 构实现，子结构的对接是精确的，不会因为子结构不具有自然部件特性而造成 计算上的误差，恰好能够符合子结构自动分解的结果。 除以上并行计算之外，论文作者攻读学位期间还在有限元结果可视化方面 进行了探索，成果已在计算机图彤和图像学报上发表，也作为学位论文一 部分内容。 全文由三个主要部分组成。 第一部分：介绍子结构自动分解 第二部分：介绍静态子结构并行计算和子空间迭代法的子结构实现 第三部分：介绍三维实体有限元模型的一种可视化方案 论文的主要创新点及其特色在于： 夺提出了新的子结构自动分解方法，能够随意控制予结构的数目和规模， 更加适应在p c 建构的松散耦合多处理机上实现并行计算。 夺设计了广义特征值问题子空间迭代法的子结构实现技术，可以得到与整 体结构完全等价的计算结果，结果精度不受子结构分解形式影响。 令提出用透明效果表现有限元模型的可视化方案，不仅可以显示三维实体 中高应力区表面分布，而且可以显示其深度，增强了可视化表现力。 夺论文内容包含子结构自动分解、静力及模态分析、可视化处理，涉及有 限元分析全过程，文中提出的各种方法都得到工程实例验证。 4 2 子结构自动分解 现代有限元分析更多是借助通用程序对结构进行网格自动生成。网格自动 生成的结果将得到一个结构完整的结点几何信息、单元拓扑信息、边界约束及 载荷信息，以及相关的材料数据等，并以特定的格式记录在有限元数据文件中。 子结构自动分解就是从数据文件中提取这些信息，然后根据一定规则自动地分 解成若干个有利于并行计算的子结构信息。 由于子结构自动分解是直接从整体有限元数据文件中提取信息的，因此在 用通用程序生成有限元网格时不必顾及子结构的分解要求，可以尽量利用通用 程序的优势，提高工作效率，给网格生成以更大的自由度。也就是说，工程师 进行有限元分析只要关心网格是否能够满足计算效果的要求，至于以后的怎样 进行子结构分解和并行计算对于他们来说是完全透明的。 一般来说，子结构的分解得恰当与否，会在很大程度上影响着基于子结构 并行计算效率。目前，有关于如何自动分解子结构的文献”“7 8 】一般都是提出 一些原则意见，具体怎样实现很少有专门报道专题。由于从本文第二部分将要 介绍的子结构静动力分析完全等价于整体结构的计算，所以计算精度与子结构 分解是无关的。因此自动分解子结构方法的选择更多是从计算效率上进行考虑。 文献t 9 1 0 ，“，1 2 1 所提出了三种区域( 子结构) 分解的算法：坐标递归二分法( r c b ： r e c u r s i v ec o o r d i n a t eb i s e c t i o n ) ，图形递归二分法( r g b ：r e c u r s i v eg r a p hb i s e c t i o n ) 和特征递归二分法( r s b ：r e c u r s i v es p e c t r a lb i s e c t i o n ) 。这三种算法都是假设计 算区域是连通的、可递归的，先按某种方法把区域( 结构) 分解成两个子区域， 然后按相同的算法作用于每个子区域，这样经过k 次后，就可得到2 个子区域。 r c b 方法【”1 l ”l 根据区域中结点的坐标来进行分解。假设v = v ，v 。v 。) 是一点集空间，v ，【j _ 1 , 2 ，月) 是空间的一个二维( x ，只) 或者三维( x ，y 。，z ，) 的 点，分解算法如下： 1 找出区域v 中结点坐标的最大值( 在五，或z 方向) 。 2 在指定的方向按照结点的坐标大小进行归类。 3 分配给每一个子结构一半的结点。 4 然后再在每一个子结构重复上述步骤。 这种算法简单明了，并且分解的速度很快。但缺点也很明显，它只考虑到 结点的几何性质，而没有考虑到空间拓扑性质，不适合处理非均匀有限网格问 题。 r g b 方法 1 3 , 1 4 , 15 l ，它考虑到了r c b 方法没有充分利用图形的连接信息的不 足，它用拓扑距离代替欧几早德f l i 离作为度量标准，算法同r c b 相似： 1 运用s p a r s p a kr c m 算法”对区域中的结点进行编号，找出两个拓扑 距离最大的结点。 2 以这两个结点作为根结点，其他的结点按与那一个根结点的距离近就 属于此根结点所形成的子区域。 3 反复使用以上算法递归下去，直到子区域大小满足条件为止。 r g b 方法的分解效果要比r c b 方法好，得到的子区域的界面自由度数要 少一些，也适合处理非均匀有限网格问题。 r s b 方法1 1 6 】1 ”，1 ”要复杂一些，它考虑了图( 即要被分解的结构) 的l a p l a c e 矩阵l ( g ) = ( 1 0 ) ，( i ，j = 1 , 2 ，”) ，其中 7 。= 节1 点f 的度数守和节村相连 ( 即与节点外目连的边数) 其z=它j0 * o 显然上( g ) = 一d + a ，其中爿是图g 的连接矩阵，而d 是结点的度数所构 成的对角矩阵。可以验证丑= 0 是矩阵( g ) 的一个特征值，显然l ( g ) 负定矩阵， 所以 = 0 也就是l ( g ) 的最大特征值。设旯：和v ：是l ( a ) 的次大特征值和对应 的特征向量，文献口”认为利用v ，可以对区域进行分解： 1 使用l a n c z o s 变换计算l ( g ) 次大特征向量( f i e d l e r 特征向量) v ，。 2 根据在f i e d l e r 特征向量中所对应的特征值的大小对结点进行分类。 3 各分一半的结点给各个子区域。 4 再对各个子区域进行递归以得到更多的子区域。 r s b 分解效果要比上述两种都好，它能很好的处理各种形状的区域。但由 于要计算l a p l a c e 矩阵工回的特征值问题，所以计算量较大。 上述三种算法各有优劣，而且最终都将结构分解成2 。子结构。但在网络并 行计算环境中工作站的数目往往不恰是2 个如果只能限定分解成2 各予结 构，那么就有可能使得部分工作站处于闲置状态。为了充分利用计算机设备， 我们在本文中提出了一种新的算法，更能适合网络并行计算环境。 2 1 子结构的分解流程 对自动网格生成的有限元数据文件进行分析，从中提取整个结构的外表面 以及内表面信息，然后从外表面的结点中确定p 个基点p 是投入计算的工作站 数目) ，再以基点作为中心按23 节所介绍的算法对整个结构实施分解。流程如图 l 。 每个子结构的分解都是从基点开始，基点选择的标准将在2 2 节中说明。分 解后的子结构信息包括所含的单元号，内部结点号，外部结点号，算法将在2 3 6 中介绍。当然为了计算需要分解子结构以后还要确定子结构之间的连接信息 细节在2 4 节中说明。 开始 1 分析整个结构的数据 确定基点 i 子结构分解 i 确定子结构连接信息 结束 图1 子结构的处理流程 2 2 基点的确定 子结构的大小不是指子结构体积( 或面积) 的大小，而是指子结构所包含的 单元多少。同样，在确定基点时需要度量结点间的距离也不是通常意义下的欧 几里德距离，而是图论意义下的拓扑距离。结点之间的拓扑距离定义为连接两 个结点的最少单元数目。以二维有限元网格为例( 图2 ) ，比如结点1 与结点6 由 单元( 1 ) 相连，所以它们之间的距离为1 。而结点1 与结点1 9 由单元( 1 ) ，( 5 ) ，( 8 ) ，( 1 1 ) 或( 1 ) ，( 4 ) ，( 7 ) ，( 1 0 ) ，( 1 1 ) 等等路径相连，取其中最小的路径，所以拓扑距离为4 。 4 3 2 ( 3 )( 6 )( 9 )( 1 2 ) 7【l1 51 9 ( 2 )( 5 )( 8 )( 1 1 ) 61 01 41 8 ( 1 )( 4 )( 7 )( 1 0 ) 5913i7 图2 二维有限元模型 2 2 1 拓扑距离计算 计算两个结点n l ，n 2 拓扑距离的算法为： 1 令d = o ： 2 若n 1 = n 2 ，返回距离d ：否则转3 。 3 构造n l 的邻点集n l i s t ，即与结点n l 同属一个单元的所有结点的集合。 4 令d = d + i 。在n l i s t 中检索是否存在1 1 2 。如存在返回距离d ，否则转5 。 5 构造n l i s t 的剩余邻点集r l i s t ，即与集合n l i s t 中的结点同属一个单元， 且未被检索过的结点集合。令n l i s t = r l i s t ，转4 。 以图2 为例，计算结点n l = l 和结点n 2 = l o 的拓扑距离。因为n l n 2 ，所 以需要构造结点l 的邻点集n l i s t = 2 ，5 ，6 ，这是d = d + i = 1 。由于n l i s t 未检索 到结点1 0 ，于是再构造n l i s t 的剩余邻点集r l i s t = 3 ，7 ，9 ，1 0 ，1 1 ) 。令n l i s t = r l i s t ， 结果在新的n l i s t 中检索到结点l o ，最终结点1 与结点1 0 的拓扑距离d = d + i = 2 。 2 2 2 基点选择 由于基点的选取直接关系到后面各个子结构的信息提取，会影响到各个子 结构的形状、大小以及边界点的容量。通常以各个子结构的大小均匀、形状规 则而边界结点少为最佳方案。所以基点的选取非常重要。我们选择基点的准则 是基点之间的最短拓扑距离达到最大。 但如果用穷举法比较整个结构的所有结点之间的拓扑距离计算量将非常庞 大。为了提高效率，基点只要在结构的边界上选择。基点选择的流程如图3 ( 以 分成三个子结构为例) 。若有必要，调整结点的过程还可以继续下去，但是实际 计算表明，仅仅经过一轮调整过程便能收到很好的效果。 提取结构的边界点信息 从中随机选两个节点n 1 ，1 1 2 在结构外表面节点中，找出第3 个节点n 3 ， 使得n 3 到n 1 和n 2 的最短拓扑距离达到极大。 调正n 1 ，使得n l 到n 2 和n 3 的最短拓扑距离达到极大。 【调正n 2 ，使得n 2 到n l 和n 3 的最短拓扑距离达到极大。 图3 提取基点的算法流程图 2 3 子结构分解 以把整个结构分成三个子结构为例，分解的过程是这样进行：根据所得的 三个基点，分别选择三个包含这些基点的单元，构成三个原始子结构。然后采 用抢占的方式，依次将剩余的单元收集为相邻的原始子结构，直至所有单元划 归到各个子结构。若遇到一个单元同时与两个或更多的原始子结构相邻，则可 以考虑划归到体积较小的子结构，这样有助于子结构的均称。如果在网络并行 环境中各个工作站的性能不一致，那么子结构的大小比例可以人为地选择，以 适应工作站的性能，使得计算过程尽可能保持同步。 具体算法如下： 1 对结构中的所有单元作标记，并令标记数组f l a g e l = o ； 2 构造单元霹f 接集e l i s t e ，e l i s t e n 与单元e 的邻接单元编号的集合。 3 根据基点信息确定原始子结构的单元号e 1 ，e 2 和e 3 ，修改单元e 1 e 2 和 e 3 的标记f l a g e 1 = l ，f l a g e 2 = 2 ，f l a g e 3 = 3 ； 4 构造e l i s t e 】的剩余单元邻接集r l i s t e ，即与集合e l i s t e 】中的单元相邻 并且未被抢占的单元，将r l i s t e 仲的元素都划归e 所形成的子结构所 有，且令f l a g e = 相应的子结构编号且e l i s t e = r l i s t e ，转5 。 5 检索标记数组f l a g e n 值，若全部不等于0 ，说明单元已全部分解完 毕。否则转4 。 至于每个子结构的内、外结点号，我们可比较任意两个子结构中的单元号， 由单元与单元之间的结点号是否相同来判断哪些结点是内结点，哪些结点是外 结点。 2 4 予结构拓扑信息提取 为了保持程序的通用性，我们对得到的每个子结构的信息采用同一种格式 的数据文件方式记录，并且，还要知道每个子结构的总接点数，内部接点数与 外部接点数，单元数以及外部接点在总体外部界面所处的位置。子结构的内部 接点与外部接点是这样定义的：若接点仅属于唯一的个子结构，则我们认为 它是内部接点，否则是外部接点。为了方便起见，我们在每个子结构的输出数 据文件中都采用内部接点排在前面，随后是外部接点。 另外，还要找出每个子结构的外部点在总体边界点中的位置。算法流程如 图4 为： 9 输入每个子结构的信息 i 获得每个子结构的外部 接点信息 i i l 提取总体边界接点信息 i 输出每个子结构的外部 点在总体边界点的位置 信息 图4 予结构拓扑信息提取框图 至此，介绍分解子结构的算法已结束了。 2 5 子结构分解实例 分解实例选用航空发动机涡轮空心叶片的有限元模型。有限元模型中的结 点数为4 7 6 9 个，单元有2 7 4 7 个，如图5 。我们用三种过程分解子结构，三种过 程的分解效果如表1 。 1 随机选取边界上的三个点作为基点，直接用2 3 节的方法分解子结构。 结果如图6 。这种分解过程得到的子结构大小相差悬殊，效果不好。 2 先随机选取边界上的三个初始点，用2 2 节的方法对初始进行调整，然 后再用2 3 节的方法分解子结构。结果如图7 。这种分解过程的效果有 显著提高，基本满足子结构分解要求。 3 对过程2 的基点再作一轮调整得到新的基点，然后再用2 3 节的方法分 解子结构。结果如图8 。比较图7 和图8 可以得知，重新调整基点取得 的效果更好，但对最终子结构的规模影响不大。 图5 叶片原型n 6 第1 种过程的分解结果 图7 第2 种过程的分解结果图8 第3 种过程的分解结果 表1 分解效果比较表 子结构序号结点总数边界结点数单元数分解效果 过 11 1 9 03 2 65 7 7 差 程 29 2 44 6 13 7 9 1 33 2 9 74 8 11 7 9 1 过 11 8 5 04 2 39 1 3 较好 程 21 8 4 73 4 49 1 6 2 31 5 6 41 9 49 1 8 过 11 5 5 41 8 39 1 6 较好 程 21 7 4 82 2 49 1 6 3 31 7 9 82 4 69 】5 3 基于子结构的并行计算 3 1 分布式网络并行计算环境的构建 所谓分布式网络并行计算坏境，即指将多台工作站或微机( 统称节点机) 通 过网络联接，再配以必要的操作系统( 如p v m ，m p i 等) 所构成虚拟并行计算机。 用户操作时并不感觉到网络的存在，如同是在操作台真实的计算机。p c 并 行群机的构建将给我国提供极为廉价的超级并行计算机，这对我国的科技进步、 经济建设和教育事业具有深远的影响。 利用实验室现有的三台微机，通过集线器联接成一个局域网，如图9 所示。 机器的配置分别为： 主机2 6 6 ：c p u - - p i l l 6 6 ，内存6 4 m b 。 从机4 5 0 ：c p u - - p i l 4 5 0 ，内存6 4 m b 。 从机1 6 6 ：c p u p e n t i 啪1 6 6 ，内存3 2 m b 。 操作系统都是w i n d o w s 9 8 。在每台机器中安装m p i 构成一个并行计算环境。 大小。 图9 小型局域网示意图 由于各台微机之间的性能相差较大，所以在分配任务时应考虑到工作量的 3 2m p i 的概况 完成参加并行计算的节点机的配置，交互以及任务之间的联系需要专门的 软件。一般有p v m ( p a r a l l e lv i r t u a lm a c h i n e ) ，m p i ( m e s s a g ep a s s i n gi n t e r f a c e ) ， e x p r e s s 。p 4 ，l i n d a 。在构建分布式网络并行计算环境时我们选择的软件平台是 m p i ，因此将以m p i 为例作详细介绍。 m p i 是为开发基于消息传递的程序而制定的标准。参与m p i 标准化工作 的主要包括来自欧美的四十多个组织的六十多个人。大部分主要的并行机制造 商、大学的研究所、政府实验室、以及工业组织都投入了m p i 的标准化工作。 其中i b m 的t jw a t s o nr e s e a r c hc e n t e r ，i n t e ln x 2 ，e x p r e s s ，n c u b ev e r t e x 1 2 一瀛 汰 譬一 一。跳 m “ 。量凇 和p a r m a c s 的工作对于m p i 的形成起了巨大的作用。另外z i p c o d e ，c h i m p ， p v m ，c h a m e l e o n 和p l c l 也有重要贡献。 创建消息传递标准的显著优点是可移植性和易用性。在低层的消息传递库 的基础上，创建分布式内存通信环境的高层函数库时，这些优点是显而易见的。 更进一步的消息通信标准的益处还在于为供应商提供了一套定义清晰的基础 库，使得他们可以更有效的加以扩充，或在某种情况下为某些功能提供硬件支 持。 m p i c h 和l a m 是通用并行消息传递接口m p i 的两个实现。m p i c h 是由 a 唱0 1 1 l l en a t i o n a ll a b o r a t o r y 和m i s s i s s i p p is t a t eu n i v e r s i t y 合作开发的。l a m 是 由o h i os t a t eu n i v e r s i t y 开发而成。 x m p i 是l a mm p i 系统的x 窗口界面，它由美国o h i os t a t eu n i v e r s i t y 开 发而成，界面采用m o t i f 工具包写成。开发x m p i 的主要目的是为m p i 系统提 供一个具有用户界面的控制台。它本身是一个普通m p i 任务，通过与t r a c e 文 件通信获得并行虚拟机中其他任务的运行、通信状态以及虚拟环境的配置信息， 并借助于x 窗口界面以图形的方式进行实时显示，这些信息使得用户能够直观 地观察、监控p v m 应用程序的运行过程，有助于应用程序的调试及通信过程 的优化处理。 m p i 是一个消息传递型并行通信的程序设计规程，它是一个消息传递库， 其库函数可以嵌在c 、c + + 和f o r t r a n 等语言中。消息传递的并行编程主要是通 过调用m p i 消息传递函数库来进行的，它实现了处理机间的的数据交换功能， 并提供了并行任务之间的同步和收发数据的接口。每个任务的发送操作必须与 一个接受操作的任务相匹配，反之亦然。这就完善地提供了同步算法和异步算 法的程序实现机制。 m p l 通信规约主要提供阻塞式( b l o c k i n g ) 和非阻塞式f n o n b l o c k i n g ) 通讯方 式，通信应答关系十分严谨，两种方式都提供点对点通讯和聚合通信 ( c o l l e c t i o n ) 。 点对点通信是最基本的通信模式，有标准( s t a n d a r d ) 、同步( s y n c h r o n o u s ) 和 预备( r e a d y ) - - 种模式。 聚合通信包含聚合同步和数据交换两方面。 聚合同步调用该操作的任务要等到组内所有成员都达到该同步点后程序 进程才继续往下执行。 数据交换在一组任务之间一起进行数据交换，包含以下四种模式： 广播( b r o a d c a s t )组内一个成员的数据发给所有成员。 数据分发( s c a t t e r )组内横向顺序数据交换转为纵向数据。 数据聚集( g a t h e r )将一个组内纵向数据交换按序传给各个组。 全数据聚集( a l l t o a l l )将所有组内数据进行纵横方向的交换。 至于m p i 程序的特点及范例在附录中将进一步说明。 3 3 并行算法与串行算法 计算机体系结构有单指令流单数据流s i s df s i n g l ei n s t r u c t i o ns t r e a m & s i n g l ed a t as t r e a m ) 、单指令流多数据流s i m d ( s i n g l ei n s t r u c t i o ns t r e a m & m u l t i p l e d a t as t r e a m ) 、多指令流多数据流m i m d ( m u l t i p l ei n s t r u c t i o n s t r e a m & m u l t i p l e d a t as t r e a m ) 、多指令流单数据流m i s d ( m u l t i p l ei n s t r u c t i o ns t r e a m & s i n g l ed a t a s t r e a m ) 等4 类。 并行计算机主要可分为两大类： 1 单指令流多数据流s 1 m d ，各处理器同一时刻执行相同的指令，但所 操作的数据不同。 2 。多指令流多数据流m i m d ，各处理器对于不同数据可以独立进行操作， m i m d 多处理器一般又分为两类：共享存储与分布式存储。 通常的程序设计所采用的串行算法体现为s i s d 。利用并行计算机进行并 行计算一般要有专门的算法。而利用网络进行并行计算，由于网络间消息传递 时间较大，一般采用粗粒度的算法。并行算法的设计有两种方法：一种是直接 从数学物理问题出发，面向并行系统研制并行算法：另一种在原有的串行算法 基础上改进。 并行程序设计的过程为 1 把一个大任务分解为若干个小任务，找出可并行部分，使用消息传递 库中的函数进行编程。 2 把这些小任务分散到不同的接点作为进程同时执行。 3 节点中的并行任务进程利用消息传递库中的函数彼此通信。 在进行并行程序设计时，要考虑到两个方面的问题：一是每次传送的数据 块的大小，因为数据块太大会造成延时，影响到整个程序的运行效率。二是总 体的数据传送量。由于计算机在内存与本地硬盘的传送速度远比在网络上传递 速度快，所以尽量减小并行计算各个节点间数据传送，这样会大大地提高并行 计算的效率。 3 4 子结构的静态凝聚 设从结构中任意取出第j 个子结构 成下列的分块形式： 粼( n i l ( ”n 心j ) n = 。j l 。jl 磴) j 它的有限单元的基本方程可经过排列 ( 1 ) 1 4 其中，上标_ ，表示子结构的编号，下标i 表示子结构的内部，b 表示子结构的边 界，缸；q 为第个子结构的内部结点位移向量，讧 为第，个子结构的边界 结点位移向量。将( 1 ) 展开成为： 愀m j l + 【础 纠) = ) 七；? “；q br l - ( j 】 “：n = 碍川) 由( 2 ) 式知，内部位移函) 可以用边界位移表示岳 ， 妇门 - b 一，】_ 1 k 一慨，】- 1 融】_ 1 船力 ( 3 ) 将( 3 ) 代入方程( 2 ) 得到子结构的边界结点力 胃1 ：“础卜【f 舶【t ：门 一1 【膏】地；q + 女；巾 t j “把“1 j ( 4 ) 引入记号： 鼢】：幽】一牌，】。1 妇】 协q ：眇 一牌， 则( 4 ) 可写成： 。 【 ：，】伽p ， = 扁) ( 5 ) ( 6 ) ( 7 ) ( 7 ) 即为子结构，的静态凝聚后的方程。 整个结构的全部子结构作上述运算，并将得到的方程( 7 ) 中各矩阵经过扩大 再一对号叠加”，组成整个凝聚后子结构的耦合平衡方程： 哦l u 。) = 哦 ( 8 ) 式中 l i d ：主刚 ( 9 ) ；l 每。 ；芝) 0 0 ) ( 9 ) 和0 0 ) 中的求和运算被理解为“对号叠加”， m 是子结构的总数。【丸】 表示缩减的整体刚度矩阵，池 为各个子结构公共边界上的结点位移向量。 3 5 基于子结构的静力问题并行算法 在2 5 节子结构自动分解的基础上，结合静态凝聚技术就可以在网络并行计 算环境中分析涡轮空心叶片在高速旋转下的应力场。考虑到我们所构建的虚拟 并行机中的节点机性能有差异，为了更合理地安排工作负荷，在分解子结构时 每个子结构的单元数作了优化配置。其中第1 个子结构单元数为5 7 7 ，结点数 位5 7 7 ，由节点机主机2 6 6 承担计算。第2 个子结构单元数为3 7 9 ，结点数位9 2 4 ， 由节点机从机1 6 6 承担计算。第3 个子结构单元数为1 7 9 1 ，结点数位3 2 9 7 ，由 节点机从机4 5 0 承担计算。算法的流程如图1 0 。节点机之间的数据传输和工程 进程如图1 1 。 下面的代码为m p i 并行程序c + + 语言伪代码： i f m y r a n k = = 0t ，我是主机 从数据文件中得出子结构数目n u m s u b s t r u c t u r e s ； 把n u m s u b s t r u c t u r e s 发送给其他从机： ) e l s e ，我是从机 接受从主机传递过来的n u m s u b s t r u c t u r e s ： 定义： c m i ，c v i子结构的刚度矩阵栽荷向量 l p i子结构装载信息以及边界结点数 n j o i n t s i子结构边界结点数 p m o d e l子结构有限元模型对象 i fm y r a n k = = 0f l 我是主机 从数据文件中读取各个子结构的信息； 分析得出i p i ，n j o i n t s i ； 把i p l ，n j o i n t s l 传送给节点机1 ： 把i p 2 ，n j o i n t s 2 传送给节点机2 ： ) e l s e 我是从机 接收来自主机的数据 提取凝聚后刚度矩阵c m i ； 提取凝聚后载荷向量c v i ； b a r r i e r ( ) ；保证各个节点机的工作保持同步 i fm y r a n k 0 我是从机 ( 把已凝聚后刚度矩阵c m i 传送回主机： 把已凝聚后刚度矩阵c v i 传送回主机： 接收已由主机进行过耦合处理的c m i ，c v i ； ) e l s e 我是主机 接收来自从机的数据c m 2 ，c v 2 和c m 3 ，c v 3 叠加本机所计算的c m l ，c v l ： 进行耦合分析： 1 6 把已进行过耦合处理的c m 2 ，c v 2 和c m 3 ，c v 3 传送回从机。 计算子结构内部位移c v i ； 计算子结构的应力s t r e s s i 输出到文件中完成计算。 将子结构有限元数据发送到各个节点机( 串行) i l各节点机生成刚度矩阵和载荷向量，并进行静态凝聚， l然后将凝聚后的刚度矩阵