一类变系数二阶椭圆交面问题有效的谱元法
An Efficient Spectral Element Method for A Class of Second Order Elliptic Interface Problems with Variable Coefficients
摘要: 本文针对圆域上一类变系数二阶椭圆交面问题提出了一种有效的谱元法。首先,根据极坐标变换公式,极条件以及傅里叶基函数展开,原问题被分解为一系列关于径向变量的相互独立的一维二阶问题,并建立了弱形式和相应的离散格式。其次,我们根据变系数的正则性,构造了基于分片高阶多项式逼近的一种谱元方法,再通过利用勒让德多项式的正交性质,构造了一组适当的基函数,使得离散变分形式中的系数矩阵在变系数为分片多项式条件下是分块对角的稀疏矩阵。最后,我们呈现了一些数值例子,通过数值结果验证了我们提出的算法是收敛的和高精度的。
Abstract: In this paper, we put forward an efficient spectral element method for a class of second order ellip-tic interface problem with variable coefficients in a circular region. Firstly, because of the polarity transformation, pole condition and fourier basis function expansion, the original problem resolve into a succession of independent one-dimensional second-order problems about radial variables. Furthermore, the weak form and relevant discrete scheme are established. Secondly, according to the regularity of variable coefficients, we construct a spectral element method based on piecewise high-order polynomial approximation, and then by taking advantage of the orthogonal property of Legendre polynomials, we construct a set of appropriate basis functions, so that the coefficient ma-trix in the discrete scheme is sparse diagonal matrix in the case that the variable coefficients are piecewise polynomials. Finally, some numerical examples are presented, and the numerical results show that the proposed algorithm is convergent and high-precision.
文章引用:何娅, 郑继会. 一类变系数二阶椭圆交面问题有效的谱元法[J]. 应用数学进展, 2023, 12(4): 1438-1445. https://doi.org/10.12677/AAM.2023.124147

1. 引言

在材料科学和流体动力学中,通常会涉及到一些分片光滑的变系数二阶椭圆交面问题的数值求解。常用的数值方法包括有限元方法 [1] [2] [3] [4] [5] 和有限差分方法 [6] [7] [8] 。

文献 [9] 应用拟合有限元方法解决椭圆界面问题,并在有限元空间的一些近似假设下,得到了能量范数估计。文献 [3] 利用具有特殊拟合网格的标准有限元对拟合网格有限元方法进行估计。文献 [10] 提出了一种无限元素方法,可以看作是一种特定的网格细化方案,用于不适合弯曲界面的椭圆界面问题。文献 [11] 利用Nitsche方法,提出了椭圆界面问题的有限元法。然而,很少有关于用谱元法去处理圆域上一类分片光滑变系数二阶椭圆交面问题的报道,因此,本文针对圆域上一类变系数二阶椭圆交面问题提出一种有效的谱元法。首先,根据极坐标变换公式,极条件以及傅里叶基函数展开,原问题被分解为一系列关于径向变量的相互独立的一维二阶资源问题,并建立了弱形式和相应的离散格式。其次,我们根据变系数的正则性,构造了基于分片高阶多项式逼近的一种谱元方法,再通过利用勒让德多项式的正交性质,构造了一组适当的基函数,使得离散变分形式中的系数矩阵在变系数为分片多项式条件下是分块对角的稀疏矩阵。最后,我们呈现了一些数值例子,通过数值结果验证了我们提出的算法是收敛的和高精度的。

本文剩余部分安排如下:在第2节,我们推导圆域内变系数二阶椭圆问题的降维格式。在第3节,我们建立相应的弱形式及其离散格式。在第4节,我们详细描述算法的实现过程。在第5节中,我们呈现了几个数值例子。

2. 降维格式

作为一个模型,我们考虑下面的分片光滑的二阶椭圆交面问题:

Δ u ^ ( x , y ) + | a ^ ( x , y ) c | u ^ ( x , y ) = f ^ ( x , y ) , ( x , y ) Ω , (1)

u ^ ( x , y ) | Ω = 0 , (2)

其中 a ^ ( x , y ) = r 是一个径向函数, Ω = { ( x , y ) : x 2 + y 2 R 2 } c ( 0 , R ) 是一个非光滑点。下面我们将推导方程(1)~(2)基于极坐标变换 x = r cos θ , y = r sin θ 的降维格式。令

u ( r , θ ) = u ^ ( x , y ) , f ( r , θ ) = f ^ ( x , y ) .

则问题(1)~(2)可以化为下面的等价形式:

1 r r ( r u ( r , θ ) r ) 1 r 2 2 u ( r , θ ) θ 2 + | r c | u ( r , θ ) = f ( r , θ ) , ( r , θ ) D , (3)

其中 D = ( 0 , R ) × ( 0 , 2 π ) 。利用傅里叶基函数展开,有

u ( r , θ ) = | m | = 0 u m ( r ) e i m θ , f ( r , θ ) = | m | = 0 f m ( r ) e i m θ . (4)

将(4)代入(3),利用傅里叶基函数的正交性可得到(1)~(2)等价的一系列一维二阶椭圆交面问题:

1 r r ( r u m r ) + m 2 r 2 u m + | r c | u m = f m , r ( 0 , R ) , (5)

u m ( R ) = 0 , ( m = 0 ) , (6)

u m ( 0 ) = u m ( R ) = 0 , ( m 0 ) . (7)

3. 弱形式和离散格式

I = ( 0 , R ) , ω = r 。定义一个通常的带权Sobolev空间:

L ω 2 ( I ) = { σ : I ω | σ | 2 d r < } ,

其中,对应的内积以及范数分别如下:

( σ , ϕ ) ω = I ω σ ϕ ¯ d r , σ ω = ( I ω | σ | 2 d r ) 1 2 .

进一步定义非一致带权Sobolev空间:

H 0 , ω , m 1 ( I ) = { u m : u m L ω 2 ( I ) , u m r L ω 2 ( I ) , u m ( R ) = 0 } , ( m = 0 ) ,

H 0 , ω , m 1 ( I ) = { u m : p u m r p L ω 2 p 1 2 ( I ) , p = 0 , 1 , u m ( 0 ) = u m ( R ) = 0 } , ( m 0 ) ,

其中,对应的内积以及范数分别如下:

( u 0 , v 0 ) 1 , ω , 0 = ( u 0 , v 0 ) ω + ( u 0 r , v 0 r ) ω , u 0 1 , ω , 0 = ( u 0 , u 0 ) 1 , ω , 0 ,

( u m , v m ) 1 , ω , m = p = 0 1 ( p u m r p , p v m r p ) ω 2 p 1 , u m 1 , ω , m = ( u m , u m ) 1 , ω , m , ( m 0 ) .

则(5)~(7)的一种弱形式为:找到 u m H 0 , ω , m 1 ( I ) ,使得

A m ( u m , v m ) = F ( v m ) , v m H 0 , ω , m 1 ( I ) , (8)

其中

A m ( u m , v m ) = I r u m v m + m 2 r u m v m + r | r c | u m v m d r , F m ( v m ) = I r f m v m d r .

P N 是次数不超过N次的多项式,定义逼近空间:

X N m = { σ m N C ( I ) : σ m N | [ 0 , c ] P N , σ m N | [ c , R ] P N } H 0 , ω , m 1 ( I ) .

则弱形式(8)相应的离散格式为:找到 u m N X N m ,使得

A m ( u m N , v m N ) = F m ( v m N ) , v m N X N m . (9)

4. 算法的有效实现

我们将在这一节详细描述算法的实现过程。首先,我们构造谱元逼近空间中的一组基函数。令

φ i ( t ) = 1 4 i + 6 ( L i ( t ) L i + 2 ( t ) ) , ( i = 0 , 1 , , N 2 ) ,

其中 t I ^ = [ 1 , 1 ] L i ( t ) 表示i次勒让德多项式。我们分别定义如下的内部基函数:

ϕ 1 , i ( x ) = { φ i ( t 1 ( x ) ) , x I 1 , 0 , x I 2 . i = 0 , 1 , , N 2 ,

ϕ 2 , i ( x ) = { 0 , x I 1 , φ i ( t 2 ( x ) ) , x I 2 , i = 0 , 1 , , N 2 ,

其中 I 1 = ( 0 , c ) I 2 = ( c , R ) t 1 ( x ) : I 1 ( 1 , 1 ) t 2 ( x ) : I 2 ( 1 , 1 ) 。则 X N m 的所有内部基函数为:

X N m = j = 1 , 2 span { ϕ j , 0 , ϕ j , 1 , , ϕ j , N 2 }

定义c处的交面基函数:

ψ 1 ( x ) = { ψ 11 = 1 2 ( 1 t 1 ( x ) ) , x I 1 , ψ 12 = 0 , x I 2 .

ψ 2 ( x ) = { ψ 21 = 1 2 ( 1 + t 1 ( x ) ) , x I 1 , ψ 22 = 1 2 ( 1 t 2 ( x ) ) , x I 2 .

则逼近空间 X N m 可表示为:

X N 0 = X N 0 span { ψ 1 ( x ) , ψ 2 ( x ) } ; X N m = X N m span { ψ 2 ( x ) } .

1) 当 m = 0 时,对离散格式(9)进行积分区间标准化:

2 c 1 1 r 1 u 0 N v 0 N d t 1 + c 2 1 1 r 1 ( c r 1 ) u 0 N v 0 N d t 1 + 2 R c 1 1 r 2 u 0 N v 0 N d t 2 + R c 2 1 1 r 2 ( r 2 c ) u 0 N v 0 N d t 2 = c 2 1 1 r 1 f 0 v 0 N d t 1 + R c 2 1 1 r 2 f 0 v 0 N d t 2 ,

其中 r 1 ( t 1 ) = 1 + t 1 2 c r 2 ( t 2 ) = 1 + t 2 2 ( R c ) + c 。将近似解 u 0 N 展开为:

u 0 N = k = 1 2 i = 0 N 2 u 0 i k ϕ k , i + k = 1 2 u 0 k ψ k .

则可得到离散格式(9)等价的矩阵形式为:

( A 1 0 0 B 1 0 B 2 0 0 A 2 0 0 B 3 0 B 1 0 T 0 c 1 c 2 B 2 0 T B 3 0 T c 3 c 4 ) U 0 = ( F 1 0 F 2 0 f 1 f 2 ) ,

其中

A 1 0 = ( a i j 01 ) , A 2 0 = ( a i j 02 ) , B 1 0 = ( b 10 0 , b 11 0 , , b 1 , N 2 0 ) T , B 2 0 = ( b 20 0 , b 21 0 , , b 2 , N 2 0 ) T , B 3 0 = ( b 30 0 , b 31 0 , , b 3 , N 2 0 ) T , F 1 0 = ( f 10 0 , f 11 0 , , f 1 , N 2 0 ) T , F 2 0 = ( f 20 0 , f 21 0 , , f 2 , N 2 0 ) T , U 0 = ( u 00 1 , u 01 1 , , u 0 , N 2 1 , u 00 2 , u 01 2 , , u 0 , N 2 2 , u 01 , u 02 ) T ,

其中

a i j 01 = 2 c 1 1 r 1 ϕ 1 , i ϕ 1 , j d t 1 + c 2 1 1 r 1 ( c r 1 ) ϕ 1 , i ϕ 1 , j d t 1 , a i j 02 = 2 R c 1 1 r 2 ϕ 2 , i ϕ 2 , j d t 2 + R c 2 1 1 r 2 ( r 2 c ) ϕ 2 , i ϕ 2 , j d t 2 , b 1 i 0 = 2 c 1 1 r 1 ψ 11 ϕ 1 , i d t 1 + c 2 1 1 r 1 ( c r 1 ) ψ 11 ϕ 1 , i d t 1 , b 2 i 0 = 2 c 1 1 r 1 ψ 21 ϕ 1 , i d t 1 + c 2 1 1 r 1 ( c r 1 ) ψ 21 ϕ 1 , i d t 1 , b 3 i 0 = 2 R c 1 1 r 2 ψ 22 ϕ 2 , i d t 2 + R c 2 1 1 r 2 ( r 2 c ) ψ 22 ϕ 2 , i d t 2 ,

c 1 = 2 c 1 1 r 1 ψ 11 ψ 11 d t 1 + c 2 1 1 r 1 ( c r 1 ) ψ 11 ψ 11 d t 1 , c 2 = 2 c 1 1 r 1 ψ 11 ψ 21 d t 1 + c 2 1 1 r 1 ( c r 1 ) ψ 11 ψ 21 d t 1 , c 3 = 2 c 1 1 r 1 ψ 21 ψ 11 d t 1 + c 2 1 1 r 1 ( c r 1 ) ψ 21 ψ 11 d t 1 , c 4 = 2 c 1 1 r 1 ψ 21 ψ 21 d t 1 + c 2 1 1 r 1 ( c r 1 ) ψ 21 ψ 21 d t 1 + 2 R c 1 1 r 2 ψ 22 ψ 22 d t 2 + R c 2 1 1 r 2 ( r 2 c ) ψ 22 ψ 22 d t 2 ,

f 1 i 0 = c 2 1 1 r 1 f 0 ϕ 1 , i d t 1 , f 2 i 0 = R c 2 1 1 r 2 f 0 ϕ 2 , i d t 2 , f 1 = c 2 1 1 r 1 f 0 ψ 11 d t 1 , f 2 = c 2 1 1 r 1 f 0 ψ 21 d t 1 + R c 2 1 1 r 2 f 0 ψ 22 d t 2 .

2) 当 m 0 时,对(9)进行积分区间标准化:

2 c 1 1 r 1 u m N v m N d t 1 + c 2 1 1 m 2 r 1 u m N v m N d t 1 + c 2 1 1 r 1 ( c r 1 ) u m N v m N d t 1 + 2 R c 1 1 r 2 u m N v m N d t 2 + R c 2 1 1 m 2 r 2 u m N v m N d t 2 + R c 2 1 1 r 2 ( r 2 c ) u m N v m N d t 2 = c 2 1 1 r 1 f m v m N d t 1 + R c 2 1 1 r 2 f m v m N d t 2 .

将近似解 u m N 展开为:

u m N = k = 1 2 i = 0 N 2 u m i k ϕ k , i + u m ψ 2 ,

则可得到离散格式(9)等价的矩阵形式为:

( A 1 m 0 B 1 m 0 A 2 m B 2 m B 1 m T B 2 m T c 5 ) U m = ( F 1 m F 2 m f 3 ) ,

其中

A 1 m = ( a i j m 1 ) , A 2 m = ( a i j m 2 ) , B 1 m = ( b 10 m , b 11 m , , b 1 , N 2 m ) T , B 2 m = ( b 20 m , b 21 m , , b 2 , N 2 m ) T , F 1 m = ( f 10 m , f 11 m , , f 1 , N 2 m ) T , F 2 m = ( f 20 m , f 21 m , , f 2 , N 2 m ) T , U m = ( u m 0 1 , u m 1 1 , , u m , N 2 1 , u m 0 2 , u m 1 2 , , u m , N 2 2 , u m ) T ,

其中

a i j m 1 = a i j 01 + c 2 1 1 m 2 r 1 ϕ 1 , i ϕ 1 , j d t 1 , a i j m 2 = a i j 02 + R c 2 1 1 m 2 r 2 ϕ 2 , i ϕ 2 , j d t 2 , b 1 i m = b 2 i 0 + c 2 1 1 m 2 r 1 ψ 21 ϕ 1 , i d t 1 , b 2 i m = b 3 i 0 + R c 2 1 1 m 2 r 2 ψ 22 ϕ 2 , i d t 2 , c 5 = c 4 + c 2 1 1 m 2 r 1 ψ 21 ψ 21 d t 1 + R c 2 1 1 m 2 r 2 ψ 22 ψ 22 d t 2 , f 1 i m = c 2 1 1 r 1 f m ϕ 1 , i d t 1 , f 2 i m = R c 2 1 1 r 2 f m ϕ 2 , i d t 2 , f 3 = c 2 1 1 r 1 f m ψ 21 d t 1 + R c 2 1 1 r 2 f m ψ 22 d t 2 .

5. 数值实验

为了验证该算法的有效性,我们将呈现一些数值实验。我们在MATLAB2017b平台上编程计算。令数值解与近似解之间的误差如下:

e ( u ^ ( x , y ) , u ^ N ( x , y ) ) = u ^ ( x , y ) u ^ N M ( x , y ) L ( Ω ) = u ( r , θ ) u N M ( r , θ ) L ( Ω ) ,

其中

u N M ( r , θ ) = | m | = 0 M u m N ( r ) e i m θ . (10)

例1我们取函数 u = ( x 2 + y 2 1 ) e ( x + y ) , R = 1 , c = 0.5 ,当取不同的M和N时,我们在表1中分别列出了数值解与近似解之间的误差结果。

Table 1. Error outcome between numerical solution and approximate solution for diverse N and M

表1. 当取不同的N和M时,数值解与近似解之间的误差结果

为了更加直观的表明算法的收敛性和高精度,在图1中依次给出了数值解与近似解的图像,在图2中依次给出了数值解与N = 30,M = 10和N = 40,M = 20时的近似解之间的误差图像。

Figure 1. Graphics of numerical solution (left) and approximate solutions (right) for N = 30 and M = 12

图1. 数值解(左)与N = 30和M = 12时的近似解(右)的绘图

Figure 2. Error Graphics between the numerical solution and the approximate solution when N = 30, M = 10 (left) and N = 40, M = 20 (right)

图2. 数值解与N = 30,M = 10 (左)和N = 40,M = 20 (右)时的近似解之间的误差绘图

表1的数据可以知道,当 N 30 , M 12 时,近似解 u m N ( x , y ) 达到了大约10−13的精度。最后,从图1图2也能得到我们提出的算法是收敛的和高精度。

参考文献

[1] Bordas, S.P., Burman, E., Larson, M.G., et al. (2018) Geometrically Unfitted Finite Element Methods and Applications: Proceedings of the UCL Workshop 2016. Springer, Berlin.
https://doi.org/10.1007/978-3-319-71431-8
[2] Bramble, J.H. and King, J.T. (1996) A Finite Element Method for Interface Problems in Domains with Smooth Boundaries and Interfaces. Advances in Computational Mathematics, 6, 109-138.
https://doi.org/10.1007/BF02127700
[3] Chen, Z. and Zou, J. (1998) Finite Element Methods and Their Convergence for Elliptic and Parabolic Interface Problems. Numerische Mathematik, 79, 175-202.
https://doi.org/10.1007/s002110050336
[4] Li, Z. (1998) The Immersed Interface Method Using a Finite Element Formulation. Applied Numerical Mathematics, 27, 253-267.
https://doi.org/10.1016/S0168-9274(98)00015-4
[5] Chou, S.H., Kwak, D.Y. and Wee, K.T. (2010) Optimal Convergence Analysis of an Immersed Interface Finite Element Method. Advances in Computational Mathematics, 33, 149-168.
https://doi.org/10.1007/s10444-009-9122-y
[6] Angelova, I.T. and Vulkov, L.G. (2007) High-Order Finite Difference Schemes for Elliptic Problems with Intersecting Interfaces. Applied Mathematics and Computation, 187, 824-843.
https://doi.org/10.1016/j.amc.2006.08.165
[7] Smith, G.D., Smith, G.D. and Smith, G.D.S. (1985) Numerical Solution of Partial Differential Equations: Finite Difference Methods. Oxford University Press, Oxford.
[8] Jarvis, D.A. and Noye, B.J. (2003) Finite Difference Solution to the Poisson Equation at an Intersection of Interfaces. ANZIAM Journal, 45, C632-C645.
https://doi.org/10.21914/anziamj.v45i0.913
[9] Babus ̌ka, I. (1970) The Finite Element Method for Elliptic Equations with Discontinuous Coefficients. Computing, 5, 207-213.
https://doi.org/10.1007/BF02248021
[10] Han, H. (1982) The Numerical Solutions of Interface Problems by Infinite Element Method. Numerische Mathematik, 39, 39-50.
https://doi.org/10.1007/BF01399310
[11] Hansbo, A. and Hansbo, P. (2002) An Unfitted Finite Element Method, Based on Nitsche’s Method, for Elliptic Interface Problems. Computer Methods in Applied Mechanics and Engineering, 191, 5537-5552.
https://doi.org/10.1016/S0045-7825(02)00524-8