1. 引言
局部先验估计是有限元方法的基本研究内容,基于局部先验估计的局部并行算法是求解偏微分方程的高效数值方法之一。2000年许进超和周爱辉结合局部亏量校正技巧的二网格离散方法,提出了有限元局部和并行算法,并用于求解二阶非对称边值问题(参见文献 [1] )。之后,许多学者发展该方法,将之成功用于Stokes方程(参见 [2] [3] ),量子特征值问题(参见 [4] ),对流扩散方程(参见 [5] ),定常不可压缩流问题(参见 [6] ),二阶椭圆特征值问题(参见 [7] [8] ),Steklov特征值问题(参见 [9] [10] )和传输特征值问题(参见 [11] )等。
流固振动Laplace模型问题来源于浸没在不可压缩流体中的一束平行管的振动,在核工程中具有相当重要的意义(参见文献 [12] [13] [14] ),其数值方法近年受到诸多学者的关注,例如Armentano等研究了该问题的hp协调有限元法及自适应算法,并给出了先验误差估计和后验误差估计(参见 [15] ),张宇等用非协调有限元方法研究了流固振动特征值的可保证下界(参见 [16] )。本文研究流固振动Laplace模型的局部先验误差估计和有限元局部计算方法,首先给出该问题的局部先验估计,基于局部先验估计建立协调有限元局部计算方案,接着对方案进行误差分析,最后给出数值实验展示局部和并行方案的效率。
关于有限元法和谱逼近的基本理论参见文献 [17] [18]。本文中字母C表示与网格尺寸无关的常数,它在不同的地方所表示的值可能不同。为方便起见用符号
表示
。
2. 预备知识
考虑下述流固振动Laplace模型(详见 [12] [14] [19] ):求
(振动频率),
(流体压力)使得
(2.1)
这里
是由流体占据的有界多边形区域,
表示
的外边界,
表示每个管与流体间的交界面,管的硬度为
,质量为m,流体视为完全不可压缩且密度为
,
是
边界上的单位外法向量。记
(参见图1)。
令
表示
上通常t阶Sobolev空间,其上范数和半范数分别为
和
,
,
装备范数
。令
,其上范数为
。记
。(2.1)的弱形式为:求
,
满足
(2.2)
其中
Figure 1. Schematic diagram of 2D area model
图1. 二维区域模型示意图
显然
是
上连续、对称、椭圆的双线性形式,
是
上非负、连续、对称的双线性形式。由文献 [19] 可知,(2.2)的解是由2K个特征对
给出的序列,其中特征值
均为正,假设特征值按递增排序:
。相应于特征值
的特征函数
使得
是一个线性无关集。
令
是
的一族三角形剖分,
是点x所在单元T的直径,
是剖分
的直径。令
表示
中位于边界
上的边组成的集合,
。令
是定义在
上的分片m次多项式空间,
。
相应于(2.2)的离散特征值问题为:求
和
满足
(2.3)
由文献 [19] 可知,离散问题(2.3)具有2K个正的特征值,假设特征值按递增排序:
,相应于特征值
的特征函数
使得
是一个线性无关集。
相应于(2.2)的源问题与相应于(2.3)的近似源问题分别如下:
求
,使得
(2.4)
求
,使得
(2.5)
记
经简单计算有
(2.6)
于是,对
有
这表明
。因此,由Neumann问题先验估计(参见文献 [20],或文献 [21] 命题4.4)可知,(2.4)的解
且满足
(2.7)
其中
,
为区域
的最大凹角。于是,对所有
,(2.2)的特征函数
。注意,对于具有至少一个交界面
的多角形区域
,必有
,因此
。
根据Lax-Milgram定理可知,(2.4)和(2.5)分别存在唯一解,于是可定义解算子
:
于是,(2.2)和(2.3)分别具有下述等价算子形式:
(2.8)
(2.9)
其中
。
由(2.7)可知
,即算子T有界。
由(2.6)有
从而
,即
也是有界算子。
令
表示所有相应于特征值
的特征函数张成的空间。
定义Ritz投影
满足
显然有
容易证明
。
利用Cea引理和Aubin-Nische技巧可得以下结论。
引理2.1 令w和
分别是(2.4)和(2.5)的解,如果
,则
(2.10)
(2.11)
(2.12)
引理2.2 算子T是紧的,且当
时
证明:可记
,
,注意到
,由(2.10),(2.11)及(2.7)有
证毕。
由谱逼近理论(参见文献 [17] 和文献 [15] 中命题3.1),有以下引理。
引理2.3 设
,对所有的
,存在正常数C和
使得当
时有
(2.13)
(2.14)
(2.15)
(2.16)
引理2.4 设
是(2.2)的特征对,则对任意
,
,有
证明:参见文献 [17] 引理9.1。
3. 局部先验误差估计
对于
,我们用符号
表示
。给定
,定义
和
分别表示
和
在G上的限制。记
,
和
。
令
,假设本文的网格和有限元空间满足以下条件(参见文献 [1] ):
(A0) 存在
,使得
(A1) 对
,存在
,使得
(A2) 对任意
,有
(A3) 对
,令
且
。对任意
,存在
使得
对
,考虑下述混合边值问题:
(3.1)
(3.1)的弱形式为:求
使得
其中
。
对(3.1),我们作下述假设。
.对任意
,存在
满足
且
.
由文献 [1] [10],我们有下述引理。
引理3.1 令
,
,
,则有
引理3.2 设(A0),(A2)及(A3)成立,且
。若
和
满足
则有
定理3.1 设
,
,(A0),(A1),(A2)和(A3)成立,则
证明:定义Ritz投影
满足
取
满足
,并取
使得在
上
且
。令
,则对任意
有
于是,由引理3.2有
因此,利用商空间等价模定理推出
证毕。
定理3.2 设
,
,(A0),(A1),(A2)和(A3)成立,则
证明:对任意
,由定理3.1有
于是
因此
证毕。
定理3.3 在定理3.2的假设下,令
是(2.3)的第j个特征对,
是(2.2)的第j个特征值,则存在
使得下列误差估计成立
证明:由(2.8)和(2.9)可推出
因此
利用定理3.2可得
证毕。
4. 局部和并行有限元方案
设
是
的一个形正规网格,网格直径
,
是包含奇点的一个子域,并设
是一个比D稍大的子域(即
)。令
是在
基础上加密得到的全局中网格,
是在
基础上通过局部加密得到的局部细网格,网格直径满足
。
方案1 (局部计算方案):
步骤1在全局粗网格
上解(2.3):求
,
使得
且
步骤2在全局中网格
上解线性边值问题:求
使得
步骤3在局部加密细网格
上解线性边值问题:求
使得
步骤4令
计算Rayleigh商
定理4.1 设
是由方案1计算所得,且
和A(0)~A(3)成立。若
(
且
),则
(4.1)
(4.2)
证明:令
满足
。由于
下面逐一估计
,
和
。为此,取
使得
。
首先,由于
及
(4.3)
取
,则可推出
(4.4)
根据方案1步骤3,有
(4.5)
令
,则对任意
有
于是
从而,由(4.3),(4.5)和引理3.2有
(4.6)
因为
联系(4.4)和(4.6),便得
(4.7)
下面利用Aubin-Nitsche共轭论证估计
。对任意
,存在
满足
设
和
满足方程
则可推出
由有限元误差估计和局部正则性假设
,有
于是,对任意
有
由此得
由(4.4)和三角不等式
有
(4.8)
利用Aubin-Nitsche共轭论证可得
将上述估计和(4.8)代入(4.7),得
联系引理2.1和2.3便得
(4.9)
类似的,由于
,我们可以得到估计式
(4.10)
下面估计
。由方案1可知
从而
因此
(4.11)
经计算有
(4.12)
取
,则可推出
于是,由(4.12)和引理3.2得
(4.13)
从而
(4.14)
由(4.4)和三角不等式
得
于是,将(4.8)代入(4.14)便得
(4.15)
因为
,由定理3.2有
则由引理2.1和2.3得
(4.16)
联系(4.16),(4.10),(4.9)和(2.10)便得(4.1)。
利用引理2.4,并注意到(2.14),便得(4.2)。证毕。
当孤立奇点个数大于1时,可以设计局部计算方案1的并行版本。
5. 数值实验
本节将给出一些数值实验来展示局部和并行计算方案的效率。数值实验是在具有2.6 GHZ CPU和8 GB RAM的DESKTOP-81QI2IO PC上,在MATLAB2015a环境中借助软件包IFEM (参见文献 [22] )进行的。我们在下述两个区域上计算(2.2)的近似特征值:
:方形腔
,中心有一个边长为
的菱形管(见图2);
:方形腔
,内部有两个边长为2的方形管(见图4)。
对
,选取
中边长为
与
的菱形管之间的区域作为局部加密区域(见图3)。对
,选取
作为局部加密区域(见图5)。我们用方案1分别计算
和
上的最小特征值的近似,然后用并行计算方案分别计算
和
上的特征值的近似值。数值结果列在表1~5中,表格中采用下述符号:
:在粗网格
上直接求解特征值问题的自由度;
:在中网格
上求解边值问题的自由度;
:在局部细网格
上求解边值问题的自由度;
:通过方案1步骤1求得的近似特征值;
,其中
是由方案1步骤2求得的解;
:由局部和并行计算方案求得的近似特征值。
文献 [15] 给出了问题(2.1)在区域
上的最小特征值的参考值
,其代数重数为2。由表1可以看到,利用局部计算方案1可以高效地计算出高精度的近似特征值。
Figure 2.
: a rhomboidal tube with edge side length
immersed in the square
图2.
:边长为
的菱形管浸没在方形腔
中
Figure 3.
: local correction domain on
图3.
:在
上的局部校正区域
Figure 4.
: two square tubes with side length 2 immersed in a square cavity
图4.
:两个边长为2的方管浸没在方形腔
中
Figure 5.
: local correction domain on
图5.
:在
上的局部校正区域
Table 1. The approximation of the smallest eigenvalue obtained by Scheme 1 on Ω 1
表1. 由方案1求得的
上最小特征值的近似值
Table 2. The approximation of the first eigenvalue on Ω 2 obtained by parallel scheme
表2. 由并行计算方案求得的
上最小特征值的近似值
Table 3. The approximation of the second eigenvalue on Ω 2 obtained by parallel scheme
表3. 由并行计算方案求得的
上第二个特征值的近似值
Table 4. The approximation of the third eigenvalue on Ω 2 obtained by parallel scheme
表4. 由并行计算方案求得的
上第三个特征值的近似值
Table 5. The approximation of the fourth eigenvalue on Ω 2 obtained by parallel scheme
表5. 由并行计算方案求得的
上第四个特征值的近似值