求解带有奇异函数泊松问题的虚拟元方法
Virtual Element Method for Solving Poisson Problems with Singular Functions
摘要: 本文运用虚拟元方法求解带有奇异函数的泊松问题。首先给出了在多边形域上的二维泊松问题,构造非协调扩展虚拟元空间,对如何处理在多边形域的拐角处出现的奇异点进行讨论,给出自由度的计算方法,引入投影算子,得到相应的虚拟元离散格式,而后对虚拟元离散形式进行误差分析,最后在多边形网格上对带角奇点的L型域问题进行数值计算,结果证实了理论分析的正确性。
Abstract: In this paper, the virtual element method is used to solve the Poisson problem with singular func-tions. Firstly, the two-dimensional Poisson problem on the polygon domain is given, the noncon-forming extended the virtual element space is constructed , how to deal with the singularity at the corner of the polygon domain is discussed, the calculation method of the degree of freedom is given, and the projection operator is introduced to obtain the corresponding virtual element discrete format, and then the error analysis of the virtual element discrete form is carried out, and finally the numerical calculation of the type domain problem with angular singularity on the polygon mesh is carried out, and the results confirm the correctness of the theoretical analysis.
文章引用:刘洋. 求解带有奇异函数泊松问题的虚拟元方法[J]. 理论数学, 2023, 13(12): 3780-3787. https://doi.org/10.12677/PM.2023.1312391

1. 引言

近年来,处理任意多边形或多面体组成的一般网格的方法受到了学者们广泛关注与研究,模拟有限差分方法 [1] ,有限元方法 [2] ,虚拟元方法 [3] [4] 等在科学和工程领域都有广泛的应用,其中,协调虚拟元方法最初由L. Beirao da Veiga等人 [3] 在2013年提出,该方法可认为是有限元方法的推广,该方法适用于多边形和多面体单元组成的网格剖分,是用稳定Galerkin离散法求解边值问题的最新进展。协调虚拟元方法已经被应用于求解许多问题,例如线弹性问题 [5] ,特征值问题 [6] ,板弯曲问题 [7] 等。相比于协调虚拟元方法,非协调虚拟元方法最初由Blance Ayuso de Dios等人 [8] 在2016年提出,非协调虚拟元方法可以降低对空间连续性的要求,非协调虚拟元方法也有广泛的应用,被应用于求解弹性问题 [9] ,中立型延迟抛物方程 [10] 等。椭圆型问题(如Poisson方程或Laplace方程)的解可能包含奇点 [11] ,这将会导致假设更规则的近似解的数值方法(比如:有限元方法)的收敛速度变差。因此,本文构造了椭圆型问题的特殊虚拟元空间,由于域的几何形状而具有奇异性。通过基于奇异函数扩展局部空间的扩展Galerkin方法,为了消除或至少减轻奇点的影响,采用的扩展Galerkin方法中,用特殊的奇异函数扩展逼近空间。本文首先给出二维泊松问题,二维泊松问题在许多工程领域的数值计算都有广泛的应用,本文利用了虚拟元空间的结构,在近似空间中隐式地插入奇异函数,它们属于问题中出现的微分算子的核。这样的奇异函数可以被插入到局部虚拟元空间中,因为这样的奇异函数被定义为有限维空间中局部问题的解,通过适当调整局部空间的边界条件,我们隐式地包含了奇异函数。最后,这样的局部空间运用非协调虚拟元方法来求解。

本文首先给出了椭圆型偏微分方程在多边形域上的模型问题,而后对如何处理在多边形域的拐角处出现的奇异点进行讨论,给出了相应的扩展非协调虚拟元空间,最后在多边形单元网格上对带角奇点的L型域问题进行了数值实验,结果证实了该方法的良好精度,并证明了该方法在H1范数和L2范数上的最佳收敛速度。

2. 非协调虚拟元方法

2.1. 预备知识

本文遵循Sobolev空间常用的表示法,给定区域 Ω 2 ,对于 k H k ( Ω ) 表示 Ω 上的Sobolev空间, ( , ) k , Ω 表示Sobolev空间标准内积, | | k , Ω k , Ω 分别表示 H k ( Ω ) 中的半范数和范数,且 k , Ω : = l = 0 k l , Ω 2 。给定 Ω R 2 是一个多边形区域 Γ = Γ N ¯ Γ D ¯ 为区域 Ω 的边界, Γ D Γ 的拓扑中的闭集, Γ D Γ N = n Ω 表示 Γ 的单位外法向量, g N Γ N 上分段光滑的函数。 { T n } n 为区域 Ω 上的一个多边形剖分序列, N V 表示 T n 的顶点的集合, ε n 表示 T n 中所有边e的集合,对于所有的 n K T n e ε n n K 表示 K 上的单位外法向量, h K 为单元的直径,h为所有单元中最大直径, ε n I , ε n B 分别表示内部边和边界边的集合。对于任意的 K T n e ε n ε n K 表示单元K上边界的集合, X K 是单元K的重心, h e 表示边e的长度, d ( A , X K ) 表示A与 X K 的距离, h Ω 表示区域 Ω 的直径,对于所有的 n ,存在一个正常数 γ ( 0 , 1 ) ,在 T n 上满足如下正则性假设:

1) 对于所有的 K 1 , K 2 T n γ h K 1 h K 2 γ 1 h K 1

2) 对于所有的 K T n ,K相对于半径大于或等于 γ h K 的球是星型的;

3) 对于所有的 e ε K h e γ h K

2.2. 模型问题

首先考虑二维泊松问题

{ Δ u = f , u Ω , n Ω u = g N , u Γ N , u = g D , u Γ D , (1)

其中, Ω 2 f L 2 ( Ω ) ,问题(1)的变分形式为:找到 u V g D : = H g D 1 ( Ω ) : = { v H 1 ( Ω ) | Γ D v = g D } ,使得

a ( u , v ) = ( f , v ) 0 , Ω + ( g N , v ) 0 , Γ N , v H 0 1 ( Ω ) . (2)

其中, a ( u , v ) = ( u , v ) ,u是解析函数与一系列和域中拐角有关的奇异项的组合。

由于在域的顶点处有奇异行为,导致问题(2)的解不是解析的,因此,为了便于说明,我们假设问题(2)的解u可以分解为:

u = u 0 + φ A (3)

其中, u 0 表示 Ω ¯ 上的解析函数, φ A 是奇异函数。接下来,我们假设 Γ N = g D = 0 ,因此问题(2)变为:找到 u V : = H 0 1 ( Ω ) ,使得

a ( u , v ) = ( f , v ) 0 , Ω , v V . (4)

2.3. 非协调扩展虚拟元空间

首先将 T n 细分为三层。第一层 T n 1 由邻近奇异顶点A的多边形组成:给定

γ ˜ > 0 , T n 1 : = { K T n | d ( A , X K ) γ ˜ h Ω } .

第二层 T n 2 是由与第一层中的单元至少共用一个边的多边形组成:

T n 2 : = { K T n \ T n 1 | K ˜ T n 1 , c a r d ( ε K ε K ˜ ) > 0 } .

第三层 T n 3 是由 T n 中剩余的单元组成:

T n 3 : = { K T n \ T n 1 T n 2 } .

类似地,将 ε n 细分为两层边。第一层 ε n 1 是由属于 T n 1 中单元边界的所有边的集合:

ε n 1 = { e ε n | K T n 1 , e K } .

第二层由剩余部分的边组成:

ε n 2 = { e ε n \ ε n 1 } .

接下来引入辅助函数和空间,对于给定的 K T n k ,首先,定义扩展函数 φ A K ( r , θ ) : = φ A ( r h K , θ ) ,在单元K上定义扩展多项式集 ˜ k ( K ) = { k ( K ) φ A K , K T n 1 k ( K ) , 。对于 e ε K ,定义扩展函数 φ A e ( r , θ ) : = φ A ( r h e , θ ) ,在边e上定义扩展多项式集 ˜ k 1 ( e ) = { k 1 ( e ) ( n e φ A | e e ) , e ε n 1 k 1 ( e ) ,

进一步,定义局部扩展虚拟元空间:

V n ( K ) : = { v n H 1 ( K ) | Δ v n k 2 ( K ) , n e v n | e ˜ k 1 ( e ) , K T n , e ε K } ,

对于所有的 v n V n ( K ) ,空间 V n ( K ) 的自由度为:

1) v n 在单元K的内部的值 1 | K | K v n m α d x α = 1 , , n k 2

2) 如果 e ε n 1 v n 在单元K边e上的值 { 1 h e e v n m ˜ α e d s , α = 0 , , k 1 , e ε K e v n ( n e φ A e ) | e d s , α = k , e ε K

如果 e ε n 2 v n 在单元K边e上的值 1 h e e v n m ˜ α e d s α = 0 , , k 1 e ε K

下面我们定义全局非协调虚拟元空间,对于 e ε n I n e ± 表示边e的单位外法向量,对于 e ε n B

n e 表示边e的单位外法向量,定义Sobolev空间

H 1 ( T n ) = K T n H 1 ( K ) = { v L 2 ( Ω ) : v | K H 1 ( K ) } ,

接下来引入跳跃算子:对于给定的 v H 1 ( T n )

v = v + n e + + v n e , e ε n I ; v = v n e , e ε n B (5)

定义子空间

H 1 , n c ( T n , k ) : = { v H 1 ( T n ) | e v n e m ˜ α e = 0 , m ˜ α e ˜ k 1 ( e ) , e ε n } .

定义全局非协调扩展虚拟元空间为:

V n : = { v n H 1 , n c ( T n , k ) | v n V n ( K ) K T n } .

2.4. 离散形式

定义 H 1 正交投影 Π ˜ : V n ( K ) ˜ k ( K ) ,对于任意 v n V n ( K ) ,满足如下正交条件,

{ a K ( v n , q ) = a K ( Π ˜ v n , q ) , K ( Π ˜ v n v n ) = 0 , q ˜ k ( K ) , (6)

由式(6)可得

a K ( v n , q ) = K v n Δ q d x + e ε K e v n q n K d s , (7)

对于所有的 e ε n ,定义 L 2 投影算子 Π ˜ 0 : V n ( K ) | e ˜ k 1 ( e ) ,对于任意 v n V n ( K ) ,满足如下正交条件,

( v n Π ˜ 0 v n , m ˜ α e ) 0 , e = 0 , m ˜ α e ˜ k 1 ( e ) , (8)

定义 L 2 投影算子 Π 0 : V n ( K ) k ( K ) ,满足 ( v n Π 0 v n , p ) 0 , K = 0 p k ( K ) Π 0 v n 可以通过内部自由度的定义来计算。

定义双线性形式 a K : H 1 ( Ω ) × H 1 ( Ω ) R ,对于任意的 u , v H 1 ( Ω ) a K ( u , v ) : = ( u , v ) 0 , K ,局部双线性形式 a ( u , v ) = K T n a K ( u , v )

定义局部离散双线性形式为

a n K ( u n , v n ) = a K ( Π ˜ u n , Π ˜ v n ) + S K ( u n Π ˜ u n , v n Π ˜ v n ) , u n , v n V n ( K ) , (9)

全局离散双线性形式 a n ( , ) : V n × V n a n ( u n , v n ) = K τ n a n K ( u n , v n ) u n , v n V n .

对所有的 n > 0 与所有 K T n ,式(9)中定义的局部离散双线性形式满足如下的一致性和稳定性。

1) 一致性:对所有 K T n 和任意的 q ˜ k ( K ) v n V n ,都有 a n K ( q , v n ) = a K ( q , v n ) 成立。

2) 稳定性:对所有 K T n 和任意的 v n V n ( K ) ,存在两个与n无关的两个正常数 α α ,有 α a K ( v n , v n ) a n K ( v n , v n ) α a K ( v n , v n )

接下来我们对右端项进行离散,在单元E定义 f n ,对于 k 2 E T n ,定义

f n , v n = E T n E f n v n d x = E T n E ( 0 f ) v n d x = E T n E f ( 0 v n ) d x ,

f h , v h 可通过自由度的定义计算。

所以,本文要解决的离散问题为,对于任意的 v n V n ,找到 u n V n ,使得

a n ( u n , v n ) = f n , v n . (10)

2.5. 误差分析

在本节中,我们对模型问题(4)的虚拟元近似(10)进行误差分析。

定义双线性形式 N n : H 1 ( Ω ) × H 1 , n c ( T n , k ) 如下

N n ( u , v ) = K T n K u n K v d s = e ε n e u v d s , (11)

其中, v 为式(5)中所定义的跳跃算子。

引理1 [12] 令u为 H 1 ( Ω ) 中的任意函数, u 0 为式(3)中所示,在 T n 上满足正则性假设(1)~(3),则在 ˜ k ( K ) 中存在分段多项式 u π ,使得

| u u π | 1 , T n c h k { ( K T n 1 | u 0 | k + 1 , K 2 ) 1 2 + ( K T n 2 T n 3 | u | k + 1 , K 2 ) 1 2 } ,

其中,c是一个依赖于k和 γ 的正常数。

引理2 [8] 令 k 1 u H s + 1 ( Ω ) , s 1 是连续问题(4)的解, v H 1 , n c ( T n , k ) N n 满足式(11)的定义式,在正则性假设成立条件下,有如下误差估计成立

N n ( u , v n ) c h min ( s , k ) u s + 1 , Ω | v n | 1 , T n ,

其中,c是一个依赖于k和 γ 的正常数.

定理1 [12] 令u与 u n 分别为(4)与(10)的解, u 0 为式(3)中所示, H 1 范数下误差估计如下

| u u n | 1 , T n c h p { ( K T n 1 u 0 k + 1 , K 2 + K T n 2 T n 3 u k + 1 , K 2 ) 1 2 + f k 1 , Ω } ,

其中,c是一个与h和u无关的常数,可能与参数 γ 以及奇异函数 φ A K 有关。

3. 数值算例

本节我们通过数值算例,验证运用非协调虚拟元方法研究带有奇异函数的泊松问题的有效性。

考虑L型域上的二维泊松方程

{ Δ u = f , u Ω , f = 0 , u Ω , n Ω u = g N , u Γ N , u = 0 , u Γ D ,

此方程的精确解为 u ( x , y ) = u ( r , θ ) = r 2 3 sin ( 2 3 θ ) ,其中, ( r , θ ) 为重入角 A = ( 0 , 0 ) 处的极坐标。

选择尺寸如图1所示的L型域

Figure 1. Schematic diagram of an L-shaped region, where A represents the vertex where the singularity occurs

图1. L型区域示意图,A表示奇点发生的顶点

在数值计算时,首先对区域 Ω 进行网格剖分,剖分区域为 Ω = ( 1 , 1 ) 2 \ ( [ 0 , 1 ) × ( 1 , 0 ] ) 图2给出了网格剖分数N是分别为53,239,834,3631,14,554,54,960时的网格剖分图。

(a) N = 53 (b) N = 239 (c) N = 834(d) N = 3631 (e) N = 14,554 (f) N = 54,960

Figure 2. Grid division diagram

图2. 网格剖分图

图3为不同网格下的误差图, e r r o r ( L 2 ) 表示 L 2 范数下的误差, e r r o r ( H 1 ) 表示 H 1 范数下的误差。可以看出随着网格剖分数的增加, H 1 L 2 范数下的误差逐渐变小。

Figure 3. The absolute errors of H 1 and L 2 norm

图3. H 1 L 2 范数绝对误差

图3给出了不同网格剖分下的 H 1 L 2 范数意义下的绝对误差。从图3中可以看出,当N变大时,数值解的误差变小,这表明了利用非协调虚拟元方法处理带有奇异函数的泊松问题的误差结果符合上述理论分析,进而证明了该方法的有效性。

4. 总结

本文运用虚拟元方法求解带有奇异函数的泊松问题。首先对该方法进行了理论分析,其次通过数值实验说明了,随着网格剖分数的增大,数值解的误差在变小,验证了理论分析结果的正确性,因此,该方法可以为其他方法提供参考,此外,扩展非协调虚拟元方法也可以用来解决更多的实际问题,例如,裂缝问题也将是我们未来的研究方向。。

参考文献

[1] Veiga, L.B., Lipnikov, K. and Manzini, G. (2014) The Mimetic Finite Difference Method for Elliptic Problems. 1st Edition, Springer-Verlag, New York.
https://doi.org/10.1007/978-3-319-02663-3
[2] Brenner, S.C. and Scott, R.L. (2008) The Mathematical Theory of Finite Element Methods. Springer-Verlag, New York.
https://doi.org/10.1007/978-0-387-75934-0
[3] Beirão da Veiga, L., Brezzi, F., Cangiani, A., et al. (2013) Basic Principles of Virtual Element Methods. Mathematical Models and Methods in Applied Sciences, 23, 199-214.
https://doi.org/10.1142/S0218202512500492
[4] Beirão da Veiga, L., Brezzi, F., Marini, L.D., et al. (2014) The Hitchhiker’s Guide to the Virtual Element Method. Mathematical Models and Methods in Applied Sciences, 24, 1541-1573.
https://doi.org/10.1142/S021820251440003X
[5] Beirão da Veiga, L., Brezzi, F. and Marini, L.D. (2013) Virtual Elements for Linear Elasticity Problems. SIAM Journal on Numerical Analysis, 51, 794-812.
https://doi.org/10.1137/120874746
[6] Lepe, F., Mora, D., Rivera, G., et al. (2021) A Virtual Element Method for the Steklov Eigenvalue Problem Allowing Small Edges. Journal of Scientific Computing, 88, Article No. 44.
https://doi.org/10.1007/s10915-021-01555-3
[7] Brezzi, F. and Donatella Marini, L. (2013) Virtual Element Methods for Plate Bending Problems. Computer Methods in Applied Mechanics and Engineering, 253, 455-462.
https://doi.org/10.1016/j.cma.2012.09.012
[8] de Dios, B.A., Lipnikov, K. and Manzini, G. (2016) The Non-conforming Virtual Element Method. ESAIM: Mathematical Modelling and Numerical Analysis, 50, 879-904.
https://doi.org/10.1051/m2an/2015090
[9] Bei Zhang, et al. (2018) The Nonconforming Virtual Element Method for Elasticity Problems. Journal of Computational Physics, 378, 394-410.
https://doi.org/10.1016/j.jcp.2018.11.004
[10] 张瑶. 中立型延迟抛物方程的虚拟元法分析[D]: [硕士学位论文]. 郑州: 河南工业大学, 2022.
[11] Grisvard, P. (1985) Elliptic Problems in Nonsmooth Domains. Pitman Pub-lishing, Boston.
[12] Artioli, E. and Mascotto, L. (2021) Enrichment of the Nonconforming Virtual Element Method with Singular Functions. Computer Methods in Applied Mechanics and Engineering, 385, Article 114024.
https://doi.org/10.1016/j.cma.2021.114024