LU 并行求解算法的研究与实现

LU 并行求解算法的研究与实现

徐 金 秀

（江南计算技术研究所 无锡 214083）

(jinxiu3186@)

摘 要 N-S 方程在CFD （computational fluid dynamics ）领域有着广泛的应用。本文介绍了仿真CFD 中N-S 方程的SSOR(symmetric successive over-relaxatio n 对称逐次超松弛)求解方法，并针对该算法数据相关性的特点给出了并行求解的不同策略及相应的分析，最后给出了一些测试结果。

关键词 并行算法、计算区域、数据分布、数据相关性、PE

1． 引言

从70年代开始世界上相继出现了多种不同类型的并行计算机，运算速度从当初的每秒200万次浮点运算达到了今天的每秒十万亿次级浮点运算。这些超级并行计算机的出现不仅为解决许多大型科学计算问题提供了硬件基础，而且对研究能充分发挥这些机器并行处理能力及解决具有重大挑战意义的科学工程计算的并行算法的发展也起到了很大的促进作用。

三维Navier-Stokes(简称N-S )方程的求解在数值天气预报、计算流体力学、石油地震勘探、航空航天等领域都有着广泛的应用，特别是由于超级计算机的出现，解决涉及数据量大、运算量大、精度要求高等问题的计算已成为可能，目前是国际上研究的热点。本文以美国NASA Ames 研究中心提出的NAS benchmark 中CFD(计算流体动力学)仿真应用LU 课题为蓝本，通过对串行程序及数据相关性的分析，采取相应的数据分配策略，从计算过程中抽取高度的并行机制，取得了较好的并行性能。

2． 数学模型 2.1 差分方程

可压缩流体与时间相关的近物面Navier-Stokes 方程可表示如下：

)1(),,,(),,,,(),(),(),()

()()(D D U U U U H U U W U U V U U T U G U F U E U ×∈+∂∂+

∂∂+∂∂+∂∂+

∂∂+∂∂=∂∂τζηξζηξζηξτζηξςη

ξτ 边界条件和初始条件分别满足下式

),,(),

,,,(),,,(),,,,(),,,(0

=∈=∂×∈=τζηξζηξζηξτςηξττςηξ当D

U U D D U U U U U B B

其中：U、E、F、G、T、V、W 和H 是5×1的向量函数；ξ、η、ζ是可计算区域的一般曲线坐标的三个方向，τ是时间变量，是有界区域，∂3ℜ∈D D 是D 的边界，}0|{

T D ≤≤=τττ。

使用一个伪时间递进步和基于有限差分近似的空间离散法可求得方程(1)的离散近似高阶耗

散项进行显式处理

时，上

述方程变为如下的线性方

程：

解。在差分过程中为了保证解的迭代收敛，要包含一个辅助高阶耗散项，当仅对辅助)2(][])()()([]})()()()()()([{*

4

444444442

22222H U h U h U h W G V F T E U S C Q B N A I n

n n n

n n n

n n n n n n τζηξετζ

ηξτζζηηξξ

τζηξ∆+∂∂+∂∂+∂∂∆−∂+∂+

∂+∂+∂+∂∆=∆∂∂+∂∂+∂∂+∂∂+∂∂+∂∂∆−

其中H ∗为修改过的向量限制函数，

上述二式左部(LHS)被称为implicit 部分，右部(RHS)被称为explicit 部分。设上述方程为如下的线性方程，

[][][]

)3(n

n

n

RHS U K =∆

其中：

[]{}

n n n n n n n S D C D Q D B D N D A D I K 222ζζηηξξτ+++++∆−=

运用中心差分对上式中的一阶偏导进行空间离散则：

)]

()()[2/1()4()]

()()[2/1()]()()[2/1(1,,1,,,1,,1,,,1,,1n k j i n k j i n n

k

j i n

k

j i n

n k j i n k j i n U C U C h C D U

B U

B h B D U A U A h A D −+−+−+−=−=−=ςςηηξξ

运用向前和向后差分对上式中的二阶偏导进行空间离散则：

)]()(2)()[/1()

5()]

()(2)()[/1()]

()(2)()[/1(1,,,,1,,22,1,,,,1,22,,1,,,,122n k j i n k j i n k j i n n k j i n k j i n k j i n n k j i n k j i n k j i n U S U S U S h S D U Q U Q U Q h Q D U N U N U N h N D −+−+−++−=+−=+−=ςζηηξξ

将(4)和(5)代入(3)式并整理得：

[]()()[]()()[]n k j i n k j i l U S h U C h A 1,,21,,/12/1−−∆−−∆−=ζζττ

[]()()[]()()[]B h B U h Q U l i j k n i j k n =−−−−−∆∆ττηη121121//,,,,

[]()()[]()()[]C h A U h N U l i j k n i j k n =−−−−−∆∆ττξξ121121//,,,,

[]()()[]

()()[]()([]

n k

j i n k

j i n k j i l U S h U Q h U N h I D ,,2

,,2

,,2/2/2/2ζ

η

ξτττ∆+∆+

∆+=

)

[]()()[]()()[]E h A U h N U l i j k n i j k n =−−++∆∆ττξξ121121//,,,,

[]()()[]()()[]F h B U h Q U l i j k n i j k n =−−++∆∆ττηη121121//,,,, []()()[]()()[]G h C U h S U l i j k n i j k n =−−++∆∆ττζζ121121//,,,,

2.3 SSOR 算法

我们将（3）式中矩阵K 写为下列形式：

K D Y Z n n n =++n 其中：

为的主对角线块。 D n K n 为的下对角线块，严格下三角。

Y n K n

Z n 为的上对角线块，严格上三角。

K n SSOR 迭代算法可写为： [][][]RHS U X n

n

=∆

其中：

()()()()()()()()

n

n n

n

n

n

n

n n n Z D I Y D

Z D D Y D X 1

1

.2.2−−++−=++−=ωωωωωωωω 2.2 算法描述

设置边界值及初始值；

、

U i j k ,,U i j k ,,0 对(i ，j ，k)，计算限制项；

∈D h H i j k ,,* SSOR 迭代：

(1) 对(i ，j ，k)，计算

∈D h []RH S i j k ,,0； (2) 循环： istep=1， etmax (a) 对(i ，j ，k)，计算

∈D h []RH S i j k n ,,；

(b) 生成并计算下面的规则、稀疏块下三角线性方程组，得到[]：

∆U 1()

[][]RHS U Y D

n n

=∆+1ω

(c) 生成并计算下面的规则、稀疏块上三角线性方程组，得到[：

]∆U n