一类Sobolev方程的各向异性非协调有限元分析

郭 城，王海红

（1.郑州师范学院 数学与统计学院，河南 郑州 450044；2.河南财经政法大学 数学与信息科学学院，河南 郑州 450002）

Sobolev方程是一类用处广泛的偏微分方程，在流体穿过裂缝岩石的渗透问题、土壤中湿气的迁移问题和不同介质的热传导问题等方面有着广泛的应用。在求解上述问题的有限元方法中，常用的有Galerkin有 限 元 方 法［1］、混 合 有 限 元 方 法［2］、H1-Galerkin混合有限元方法［3］和连续时空有限元方法［4］等，但以上研究用到的都是协调有限元方法，或要求网格剖分满足正则性假设［5］。在解决窄边区域上的偏微分方程问题时，为了减少计算量，常常采用各向异性剖分，这样可通过较少的自由度得到同样的估计结果。与协调元相比，对于自由度定义在单元的边上和单元自身上的非协调元来说，每个未知量只涉及两个单元，因此在信息传递上是廉价的，做并行计算也是比较容易的。石东洋等［6］利用积分恒等式和插值后处理技术对各向异性网格上Sobolev方程的解做了Carey元的高精度分析，但没有通过数值算例进行验证。在本文中，我们在各向异性网格上研究了非协调Carey有限元［7–8］的收敛性，在半离散格式下证明了离散问题的解的存在唯一性，给出了零模和能量模意义下的最优误差估计，通过数值算例验证了结论的正确性。

对于方程（1），我们做如下假设： Ω⊂R2是有界凸区域，a ( x, y )和 b(x,y)是具有二阶连续导数的函数； 0 < a1< a( x, y ) < a2,|b( x, y) |< a2,且 a1和 a2为正常数。

1 半离散格式下的误差估计

为了得到误差估计结果，引入经典的与时间有关的Sobolev空间［10］和相应的Hilbert空间 H，记

和

假设区域 Ω上三角形剖分族 Jh中单元 K的3条边为 li( i = 1,2,3),3个顶点坐标为 Ai( i =1,2,3),相应的顶点面积坐标为 λi，其面积为 S。

定义单元K上的形函数

在本文中，用到的有限元空间［7］为

其中,FK为标准单元K^ 到一般单元K的仿射变换,且有

由{λ1,λ2,λ3,λ1λ2+λ2λ3+λ3λ1}生 成 的 形 函 数 空间为P2,K。假设单元满足最大角条件和坐标系条件,但不需满足正则性假设或拟一致假设。由文献［7］可知，单元K具有各向异性特征。

为了便于讨论，分别定义零模和能量模

引理1［8］：对任意的有

其中(x,y)∈Ω。

式（2）的半离散格式表示如下：对于任意的vh∈Vh，求uh∈Vh，使得

其中(x,y)∈Ω。

定理1：设uh为式（3）解，则uh具有稳定性，即有成立。

证明：在式（3）中取vh=uht，由Cauchy-Scwartz不等式可得

对上式两端关于时间从0到t进行积分,由b(x,y)的 有界性 可得

定理2：设u和uh分别是式（2）和式（3）的解，则

和

成立。

证明：由引理1的证明过程可知

于是，有

由式（4）、式（5）和式（6）即可得结论。

该结论对矩形网格上的Wilson元［6］、任意四边形类Wilson元［11］以及类Carey元［12］也成立。

2 数值实验

为了验证本文结论的正确性，取Ω= [0,1]×[0,1],a(x,y) =b(x,y) = 1,f= (1 +42)etsin(x)sin(y),T=1。可 验证方程的真解u=etsin(x) sin(y)。

为了便于比较，我们采用了如图1和图2所示的两种网格对区域Ω进行剖分。其中，网格1是由沿x轴、y轴方向等分的平行线和对角线形成的平行线型三角形网格，网格2是由沿x轴、y轴方向的比为1∶8的平行线和对角线形成的平行线型三角形网格。

图1 网格1

表1给出了在网格1剖分下，时间分别为0.1、0.5和1.0时的两种类型模的误差估计。表2给出了网格1上两种模的收敛阶数。

表1 在网格1上不同变化时间的计算结果

表2 在网格1上的收敛阶

表3给出了在网格2剖分下，时间分别为0.1、0.5以及1.0时两种类型模的误差估计。表4给出了网格2上两种模的收敛阶数。

表3 在网格2上不同变化时间的计算结果

表4 在网格2上的收敛阶

从表1～表4可以看出：误差结果随着时间的变化是稳定的，在两种网格下的误差收敛阶数总体上是相同的；当h→ 0时， ||u−uh||和 |u−uh|h分别收敛于最优阶O(h2)和O(h)。由此可知，数值计算的结果与理论分析是完全一致的，这说明在传统有限元方法中，对网格剖分的正则性要求不是必要条件。