水下大型结构振动声辐射预报技术引论：并行与串行算法对比

非全局迭代法、重叠型区域分解算法的研究目标都是求解耦合方程式（2-18），而正是由于式（2-18）的复杂性质，导致了它们在收敛速度、预条件方面的困难。与它们不同的是，ODDM不再将结构有限元方程与声场边界元方程看作一个整体，而是满足位移、声压平衡条件的结构和声场两个子系统，由于两个子系统间存在着位移、声压等共同的未知量，因此称之为重叠型区域分解算法。

对结构有限元网格节点进行分类，得到外部耦合面上的节点集合和内部节点集合，分别记为O和I，用下标F和B代表有限元和边界元，得到新的声振预报方程如式（4-14）所示，位移和声压平衡条件如式（4-15）所示。

在（4-14-1）中，利用已知界面声压PF和激励力fI、fO求得位移分布向量uI，B和uO，B的问题被称为Neumann问题，利用已知界面位移T.uO，B和激励力fI、fO求解位移分布向量uI，B和pF的问题被称为Dirichlet问题；而对于声场的边界元方程（4-14-2），振速边界条件被称为Dirichlet边界条件，利用已知位移uO，B求解边界声压的问题被称为Dirichlet问题，声压边界条件又被称为Neumann边界条件，而利用已知声压pB求解边界法向振速的问题被称为Neumann问题。

分别求解（4-14）中的两个子问题，并通过两个子问题之间存在重叠关系的未知量，也就是界面上的未知量的平衡性，使得两个子问题的解向整个耦合系统的解逼近，这样的思想被称为重叠式区域分解。根据求解的子问题类型，重叠式区域分解的实现方法可以分为Dirichlet- Dirichlet方法、Dirichlet-Neumann和Neumann-Neumann方法三类。然而在本文讨论的方程形式下，界面位移向量uO，B的长度要大于边界元方程中的自由度数，我们显然无法利用已知的pB通过（4-14-2）求得完整的界面位移分布向量uO，B，要求解声场域的Neumann问题，就必须将界面上的位移表达方式由笛卡尔坐标系变换为基于法向的坐标系，加大了算法的复杂程度，不利于实现。

进一步注意到，有限元方程的Neumann问题和边界元方程的Dirichlet问题，其求解目标具有统一性，可以为平衡条件在算法中的实现提供便利，因此首先考虑一种Dirichlet-Neumann方法，步骤如算法4.1所示。

算法4.1 并行Dirichlet-Neumann算法

1）设定初始边界位移猜测向量uB0，初始边界声压猜测向量pF0，输入收敛判定值e，设定迭代计数k=0；

2）将uBk带入式（4-14-2），求解声场Neumann问题，得到边界声压向量pBk；

3）将pFk带入式（4-14-1），求解有限元方程Dirichlet问题，得到内部和界面位移分布向量uF，Ik和uF，Ok；

4）为了确保或加速迭代收敛，取松弛因子a，利用式（4-16）计算下一个迭代步的输入条件uBk+1；(www.xing528.com)

5）取pFk=pBk，根据式（4-17）计算残值εk，若εk≤e，则收敛得到满足，终止循环并返回结果向量，若未收敛，则使循环计数k=k+1，跳转至步骤2）继续循环。

算法4.1被称为并行Dirichlet-Neumann算法，是因为占据了计算过程中绝大部分计算量的步骤2和步骤3之间并不存在依赖关系，因此具备较好的可并行性。然而在现实中，各个子域问题的规模比例多变，有限元、边界元法的求解计算既受到其自身的算法复杂度影响，也受到各子域所采用的求解算法的计算复杂度影响，因此在算法4.1中，步骤2和步骤3之间的并行很难实现负载均衡。

从经典迭代算法的知识出发，可以很容易认识到，算法4.1更新改进解的方式与Jacobi迭代法十分相似，都是由上一步的整体解推算出下一步的整体改进解。而Gauss-Seidel迭代算法在每一个迭代步中，利用前j-1个已经更新过的分量替代上一步解的相应分量来更新第j个分量，在使得更新过程中不同的分量出现依赖关系、降低了算法的可并行性的同时也有效地提高了迭代收敛的速度。如上文讨论，不同子域的求解很难达到负载均衡，因此可以通过牺牲并行性而提高收敛速度来提升算法的整体计算效率。根据Gauss-Seidel迭代的思想，我们可以在算法4.1的基础上得到串行Dirichlet-Neumann算法，其具体步骤如算法4.2所示。

算法4.2 串行Dirichlet-Neumann方法

1）设定初始边界位移猜测向量uB0，输入收敛判定值e，设定迭代计数k=0；

2）求解式（4-14-2）的Neumann问题，得到边界声压向量pBk；

3）令pFk=pBk，带入式（4-14-1），求解有限元方程Dirichlet问题，得到内部和界面位移分布向量uF，Ik和uF，Ok；

4）为了确保或加速迭代收敛，取松弛因子g，利用式（4-18）计算下一个迭代步的输入条件uBk+1；

5）根据式（4-17）计算残值εk，若εk≤ε，则收敛得到满足，终止循环并返回结果向量，若未收敛，则使循环计数k=k+1，跳转至步骤2）继续循环。