Research article Special Issues

Conformable finite element method for conformable fractional partial differential equations

  • The finite element (FE) method is a widely used numerical technique for approximating solutions to various problems in different fields such as thermal diffusion, mechanics of continuous media, electromagnetism and multi-physics problems. Recently, there has been growing interest among researchers in the application of fractional derivatives. In this paper, we present a generalization of the FE method known as the conformable finite element method, which is specifically designed to solve conformable fractional partial differential equations (CF-PDE). We introduce the basis functions that are used to approximate the solution of CF-PDE and provide error estimation techniques. Furthermore, we provide an illustrative example to demonstrate the effectiveness of the proposed method. This work serves as a starting point for tackling more complex problems involving fractional derivatives.

    Citation: Lakhlifa Sadek, Tania A Lazǎr, Ishak Hashim. Conformable finite element method for conformable fractional partial differential equations[J]. AIMS Mathematics, 2023, 8(12): 28858-28877. doi: 10.3934/math.20231479

    Related Papers:

    [1] Hansaem Oh, Gwanghyun Jo . Physics-informed neural network for the heat equation under imperfect contact conditions and its error analysis. AIMS Mathematics, 2025, 10(4): 7920-7940. doi: 10.3934/math.2025364
    [2] Pengfei Zhu, Kai Liu . Numerical investigation of convergence in the L norm for modified SGFEM applied to elliptic interface problems. AIMS Mathematics, 2024, 9(11): 31252-31273. doi: 10.3934/math.20241507
    [3] Hyun Geun Lee, Youngjin Hwang, Yunjae Nam, Sangkwon Kim, Junseok Kim . Benchmark problems for physics-informed neural networks: The Allen–Cahn equation. AIMS Mathematics, 2025, 10(3): 7319-7338. doi: 10.3934/math.2025335
    [4] Yu Xu, Youjun Deng, Dong Wei . Numerical solution of forward and inverse problems of heat conduction in multi-layered media. AIMS Mathematics, 2025, 10(3): 6144-6167. doi: 10.3934/math.2025280
    [5] Ruitao Liu, Wanshan Li . A new ADI-IIM scheme for solving two-dimensional wave equation with discontinuous coefficients. AIMS Mathematics, 2024, 9(11): 31180-31197. doi: 10.3934/math.20241503
    [6] Ahmed G. Salem, Turki D. Alharbi, Abdulaziz H. Alharbi, Anwar Ali Aldhafeeri . Impact of a spherical interface on a concentrical spherical droplet. AIMS Mathematics, 2024, 9(10): 28400-28420. doi: 10.3934/math.20241378
    [7] Dario Bambusi, Simone Paleari . A couple of BO equations as a normal form for the interface problem. AIMS Mathematics, 2024, 9(8): 23012-23026. doi: 10.3934/math.20241118
    [8] Min Wu, Jiali Zhou, Chaoyue Guan, Jing Niu . A numerical method using Legendre polynomials for solving two-point interface problems. AIMS Mathematics, 2025, 10(4): 7891-7905. doi: 10.3934/math.2025362
    [9] Xin Zhao, Xin Liu, Jian Li . Convergence analysis and error estimate of finite element method of a nonlinear fluid-structure interaction problem. AIMS Mathematics, 2020, 5(5): 5240-5260. doi: 10.3934/math.2020337
    [10] Maxime Krier, Julia Orlik . Solvability of a fluid-structure interaction problem with semigroup theory. AIMS Mathematics, 2023, 8(12): 29490-29516. doi: 10.3934/math.20231510
  • The finite element (FE) method is a widely used numerical technique for approximating solutions to various problems in different fields such as thermal diffusion, mechanics of continuous media, electromagnetism and multi-physics problems. Recently, there has been growing interest among researchers in the application of fractional derivatives. In this paper, we present a generalization of the FE method known as the conformable finite element method, which is specifically designed to solve conformable fractional partial differential equations (CF-PDE). We introduce the basis functions that are used to approximate the solution of CF-PDE and provide error estimation techniques. Furthermore, we provide an illustrative example to demonstrate the effectiveness of the proposed method. This work serves as a starting point for tackling more complex problems involving fractional derivatives.



    A class of problems that give rise of singular behavior are elliptic interface problems with imperfect contact [1], and it is featured by the implicit jump condition imposed on the imperfect contact interface, and the jumping quantity of the unknown is related to the flux across the interface. The elliptic interface problems with imperfect contact have been applied to model the Stefan problem of the solidification process and crystal growth, composite materials, multi-phase flows [2,3] and the problem of temperature discontinuity between the gas and cooling solid surface [4]. More examples include the heat conduction between materials of the different heat capacities and conductivities and interface diffusion processes [5,6], the temperature discontinuity between a gas and cooling solid surface [4], the conjugate heat transfer problem in thermodynamic processes between materials that are thermally coupled through non-adiabatic contacts [7], etc.

    We consider the one-dimensional (1D) elliptic interface problem with imperfect contact

    {(β(x)ux)x+c(x)u=f(x),x(0,α)(α,1),u(0)=u0,u(1)=u1, (1.1)

    together with the following implicit jump conditions across the interface α

    {[u]=u+u=λβ+u+n,[βun]=0, (1.2)

    where 0<α<1, n is a unit normal to the interface pointing from Ω to Ω+. g± denotes the right and left limits of the function g at the point α, and [g]=g+g. Without loss of generality, we assume that Ω is computational domain, and the interface α separates Ω into two sub-domains Ω and Ω+.

    The solution of the elliptic interface problems is often discontinuous due to discontinuous coefficients or singular sources across the interface. To recover the numerical accuracy near the interface, a variety of methods have been developed via enforcing the jump conditions (1.2) and (1.1) into numerical discretization such that accurate and robust numerical algorithms can be designed. Finite difference methods (FDMs) constitute a commonly used approach for elliptic interface problems. Since the publication of the original immersed interface method (IIM) [8], they have been applied to various problems, such as the Stokes flow with elastic boundaries or surface tension [9], incompressible flow based on the Navier-Stokes equations with singular source terms [10], and nonlinear problems in magneto-rheological fluids [11]. Many further developments and analysis in various aspects of the FDM for elliptic interface problems were carried out to improve the accuracy, stability or efficiency [12,13,14,15,16,17,18,19]. Many other elegant methods have been proposed in the past decade, including the ghost fluid method [20,21], finite volume method [7,22,23,24], the matched interface and boundary (MIB) method [25,26,27], etc. Finite element methods (FEMs) constitute another common practice to resolve the elliptic interface problems. With the help of body-fitted unstructured meshes, FEMs [28,29,30] have been developed for handling elliptic interfaces and irregular geometries. For interfaces with complex topologies, the construction of high quality body-fitted meshes could be difficult and time-consuming. This motivates the development of immersed FEMs based on non-body-fitted structured meshes [31,32,33,34,35,36,37,38], etc. Nevertheless, most interface schemes in the literatures are designed to be of second-order accuracy.

    The strategy for generating higher-order difference schemes can be roughly divided into two categories. (1) The first category is expanding the stencil and including more points in the schemes [15,16,25,26,27,39,40,41,42]. Advantages of this approaches is achieving high-order accuracy and using lower order jump conditions only. The obvious disadvantages of such approach is creating large matrix bandwidths and complicating the numerical treatment near the boundaries. Gibou and Fedkiw [16] introduced an O(h4) accurate finite difference discretization for the Laplace and heat equations on an irregular domain. Previous scholars [25,26,27] presented high-order MIB methods for solving elliptic equations with discontinuous coefficients and singular sources on Cartesian grids. This type of method is based on the use of fictitious points to achieve high-order accuracy. Zhong [40] presented a high-order IIM with general jump conditions by employing fictitious points to achieve high-order accuracy and using lower-order jump conditions only. Similar to the MIB approach, a wide stencil with more points are involved at irregular points. Feng and Zhao [42] introduced a new Cartesian grid FDM based on the fourth-order accurate MIB method. Colnago et al. [15] presented a high-order IIM for solving Poisson equations with discontinuous coefficients on Cartesian grids, it combines the FDM and ghost node strategy and requires only the ordinary jumps of the function. (2) The second category involves using the differential equation and the interface relations as additional identities [13,39,41]. The advantage of this approach is that the scheme is compact and minimal stencil is used. The disadvantage is that to maintain high-order accuracy, it requires the knowledge of jump conditions of high-order derivatives. Li and Ito [13] constructed a fourth-order IIM by computing high-order jump conditions involving mixed derivatives. Linnick and Fasel [39] presented a high-order IIM for simulating unsteady incompressible flow in an irregular domain. Instead of using analytical jump conditions, they compute the jump conditions for higher derivatives numerically. Angelova and Vulkov [41] presented high-order compact FDM for elliptic equations with intersecting interfaces by using the differential equation and the jump (interface) relations as additional identities, which can be differentiated to eliminate higher-order local truncation errors. For FEMs, the higher-order convergence crucially depends on how well the interface is resolved by the triangular mesh. In practice, subparametric, isoparametric or superparametric elements are usually employed to secure the optimal order of p+1 in the L2 norm for a polynomial order of p, for both continuous [33,43] and discontinuous Galerkin [34,35,44,45,46] FEMs. It is worth pointing out that most of the numerical methods in literatures are aimed at the elliptical interface problem with homogeneous jump conditions ([u]=0,[βun]=0) or nonhomogeneous jump conditions where the jumps of the temperature as well as the conductive heat flux along the interface are known explicitly (say [u]=g1, [βu/n]=g2, with known g1 and g2). When the elliptic interface problem has an implicit jump condition as in (1.2), there are comparatively few numerical methods for solving such problems.

    In this paper, we consider the 1D elliptic interface problems with imperfect contact, where the implicit jump condition is a setting for the interface and the jump of the primary variable is proportional to the normal flux across the interface. A class of high-order FDMs is constructed for the 1D elliptic interface for both the body-fitted and non-body-fitted mesh. For each case, the second-, third- and fourth-order approximations of the implicit jump condition are provided by using the original variables on both sides of the interface and a set of jump conditions and its high order derivatives. Numerical examples are presented to verify the performance of the scheme. The numerical results show that the presented schemes can reach the theoretical accuracy for solving the elliptic interface problems with imperfect contact.

    The rest of this paper is organized as follows. In Section 2, we formulate the scheme for 1D elliptic equations with imperfect interfaces. Section 3 constructs the approximation of the normal derivative for the interface. Section 4 computes the primary variables u+ and u for the interface. Numerical examples are provided in Section 5 to demonstrate the accuracy and stability of the presented scheme. A concluding remark is given in Section 6.

    We assume that β(x) in (1.1) is a piecewise smooth function with a jump at the interface α, and that β(x) has upper and lower boundries:

    0<βminβ(x)βmax, (2.1)

    where βmin and βmax are two constants. The source term f(x) is piecewise smooth.

    Introduce a uniform grid xi=ih,i=0,1,,n with h=1n. The fourth-order compact difference scheme is constructed and three-point computational stencil for each point is adopted. For node i, the three points in computational stencil are {xi1,xi,xi+1}. Without loss of generality, assume the interface α is located in a grid interval, xkαxk+1. Because of the interface is located in the stencil, the points k and k+1 are named as irregular points, and ik,k+1 are name as regular points.

    At a regular point i, i.e., ik,k+1, for convenience, we briefly give the construction process of the fourth order compact difference scheme, for more detailed information one can refer to Ref.[47].

    βδ2xuiβxδxui+ciui+τ=fi, (2.2)

    where

    δ2xui=ui+12ui+ui1h2,δxui=ui+1ui12h,

    and

    τ=βxh23!3ux3|i+βh2124ux4|i+O(h4).

    Dropping the term τ in (2.2), it is the standard second-order central difference scheme. In order to get fourth-order scheme, we have to handle the high-order derivative terms 3ux3 and 4ux4. According to the governing equation (1.1), there is

    uxx=1β(βxux+cuf). (2.3)

    Taking the first and second derivatives to (2.3) respectively, we have

    uxxx=ϕ1uxx+ϕ2ux+ϕ3u+ϕ4f+ϕ5fx, (2.4)

    and

    uxxxx=ψ1uxx+ψ2ux+ψ3u+ψ4, (2.5)

    where

    ϕ1=βxβ,ϕ2=β2xβ2βxxβ+cβ,ϕ3=cβxβ2,ϕ4=βxβ2,ϕ5=1β,ψ1=ϕ21+ϕ1,x+ϕ2,ψ2=ϕ1ϕ2+ϕ2,x+ϕ3,ψ3=ϕ1ϕ3+ϕ3,x,ψ4=(ϕ1ϕ4+ϕ4,x)f+(ϕ1ϕ5+ϕ4+ϕ5,x)fx+ϕ5fxx,ϕ1,x=β2xββxxβ2,ϕ2,x=2βxβββxxβ2β2ββxxxβxβxxβ2cβxβ2,ϕ3,x=c(2β2xββxx)β3,ϕ4,x=ββxx2β2xβ3,ϕ5,x=βxβ2.

    Substituting (2.4) and (2.5) into (2.2) and approximating the first- and second-derivatives with the central difference, we can get the fourth-order compact scheme:

    Aiui+1+Biui+Ciui1=Fi+O(h4), (2.6)

    where

    Ai=β+ϕ1βxh26+ψ1βh212h2+βx+ϕ2βxh26+ψ2βh2122h,Bi=2β+ϕ1βxh26+ψ1βh212h2+(c+ϕ3βxh26+ψ3βh212),Ci=β+ϕ1βxh26+ψ1βh212h2βx+ϕ2βxh26+ψ2βh2122h,Fi=fiβxh26(ϕ4fi+ϕ5fx,i)βh212ψ4.

    When the interface is in contact with one of the mesh nodes, it is referred to as the body-fitted case. Without loss of generality, we assume that the interface is located the mesh node k, as shown in Figure 1.

    Figure 1.  Interface fitted with body.

    In this case, the irregular nodes are k1 and k+1, and the compact stencil of these two irregular nodes are shown in Figure {2}. Then, the high order compact scheme at these two irregular nodes can be given as

    ˉAk1uk2+ˉBk1uk1+ˉCk1u=ˉFk1, (2.7)

    and

    ˉAk+1u++ˉBk+1uk+1+ˉCk+1uk+2=ˉFk+1. (2.8)

    Apparently, there are two auxiliary qualities u and u+ in the schemes. To apply the above scheme, it is necessary to approximate the u and u+ with adequate accuracy. We prepare to use the interface jump connecting conditions and establish the linear systems with u+ and u as variables. The difficulty lies in the high-order discretization of the first-order derivatives u+x and ux on both sides of the interface.

    {[u]=λβ+u+x,[βun]=0. (2.9)
    Figure 2.  Irregular point stencil for interface-fitted mesh.

    Expand the uk+1 and uk1 at the interface α by using the Taylor series, it reads

    uk+1=u++hux|++h22!2ux2|++h33!3ux3|++h44!4ux4|++h55!5ux5|++O(h6), (2.10)

    and

    uk1=uhux|+h22!2ux2|h33!3ux3|+h44!4ux4|h55!5ux5|+O(h6). (2.11)

    Rewritting the above two formulas, we can get the representations of ux|+ and ux| as follows

    ux|+=uk+1u+hh2!2ux2|+h23!3ux3|+h34!4ux4|+h45!5ux5|++O(h5), (2.12)
    ux|=uuk1h+h2!2ux2|h23!3ux3|+h34!4ux4|h45!5ux5|+O(h5). (2.13)

    According to the governing equation, we have

    u+xx=1β+(β+xu+x+c+u+f+),uxx=1β(βxux+cuf). (2.14)

    Keeping the first two terms on the right-hand side of (2.12) and substituting (2.14) into (2.12), the second-order accurate approximate scheme of ux|+ can be given as

    ux|+=ρr,1u++ρr,2uk+1+ρr,3+O(h2), (2.15)

    where

    ρr,1=(1h+h2c+β+)/Dr,ρr,2=1h/Dr,ρr,3=h2f+β+/Dr,Dr=(1h2β+xβ+).

    To obtain the higher-order approximate of ux|+, it is necessary to handle the third- and fourth-order derivatives in (2.12). From (2.4) and (2.5), the expressions of u+xxx and u+xxxx can be given as

    u+xxx=Φ+1u+x+Φ+2u++Φ+3, (2.16)

    and

    u+xxxx=Ψ+1u+x+Ψ+2u++Ψ+3, (2.17)

    where

    Φ+1=β+xβ+ϕ+1+ϕ+2,Φ+2=c+β+ϕ+1+ϕ+3,Φ+3=(c+β+ϕ+1+ϕ+4)f++ϕ+5f+x,Ψ+1=β+xβ+ψ+1+ψ+2,Ψ+2=c+β+ψ+1+ψ+3,Ψ+3=f+β+ψ+1+ψ+4.

    Keeping the first three terms on the right-hand side of (2.12), we establish the third-order accurate approximate of ux|+ as follows:

    ux|+=ˉρr,1u++ˉρr,2uk+1+ˉρr,3+O(h3), (2.18)

    where

    ˉρr,1=(1h+h2c+β++h26Φ+2)/ˉDr,ˉρr,2=1h/ˉDr,ˉρr,3=(h2f+β+h26Φ+3)/ˉDr,ˉDr=1h2β+xβ++h26Φ+1.

    Further, the fourth-order approximate of ux|+ can be given as

    ux|+=ˉˉρr,1u++ˉˉρr,2uk+1+ˉˉρr,3+O(h4), (2.19)

    where

    ˉˉρr,1=(1h+h2c+β++h26Φ+2+h324Ψ+2)/ˉˉDr,ˉˉρr,2=1h/ˉˉDr,ˉˉρr,3=(h2f+β++h26Φ+3+h324Ψ+3)/ˉˉDr,ˉˉDr=1h2β+xβ+h26Φ+1+h324Ψ+1.

    Similarly, the second-, third- and fourth-order approximation of ux| can be given as

    ux|=ρl,1u+ρl,2uk1+ρl,3+O(h2), (2.20)
    ux|=ˉρl,1u+ˉρl,2uk1+ˉρl,3+O(h3), (2.21)
    ux|=ˉˉρl,1u+ˉˉρl,2uk1+ˉˉρl,3+O(h4), (2.22)

    where

    ρl,1=(1h+h2cβ)/Dl,ρl,2=1h/Dl,ρl,3=h2fβ/Dl,Dl=1+h2βxβ,ˉρl,1=(1h+h2cβh26Φ2)/ˉDl,ˉρl,2=1h/ˉDl,ˉρl,3=(h2fβh26Φ3)/ˉDl,ˉDl=1+h2βxβ+h26Φ1,ˉˉρl,1=(1h+h2cβh26Φ2+h324Ψ2)/ˉˉDl,ˉˉρl,2=1h/ˉˉDl,ˉˉρl,3=(h2fβ+h26Φ3+h324Ψ3)/ˉˉDl,ˉˉDl=1+h2βxβh26Φ1+h324Ψ1.

    When the interface is located within one of the computational mesh, it appears as shown in Figure 3. The irregular nodes are k and k+1, and the construction of the higher order scheme is more complicated. Owing to the interface arbitrarily cutting through the computational mesh, the computational mesh of these irregular nodes is nonuniform, as shown in Figure 4.

    Figure 3.  Interface unfitted with a body.
    Figure 4.  Nonuniform mesh.

    We assume the function u(x) is smooth enough, expand the u(xi+1),u(xi1) at point xi by Taylor series

    ui+1=ui+hfux|i+h2f2!2ux2|i+h3f3!3ux3|i+h4f4!4ux4|i+h5f5!5ux5|i+O(h6f), (2.23)

    and

    ui1=uihbux|i+h2b2!2ux2|ih3b3!3ux3|i+h4b4!4ux4|ih5b5!5ux5|i+O(h6b). (2.24)

    From (2.23) and (2.24), we can obtain

    2ux2|i=2hbui+1+hfui1(hf+hb)uihfhb(hf+hb)hfhb33ux3|i112h3f+h3bhf+hb4ux4|i160(h2f+h2b)(hf+hb)5ux5|i+O(h5f+h5bhf+hb), (2.25)

    and

    ux|i=h2bui+1h2fui1+(h2fh2b)uih2bhf+h2fhbhfhb63ux3|ihfhb(hfhb)244ux4|i+O((h3f+h3b)hfhbhf+hb). (2.26)

    Define the difference operators

    δxui=h2bui+1h2fui1+(h2fh2b)uih2bhf+h2fhb, (2.27)
    δ2xui=2hbui+1+hfui1(hf+hb)uihfhb(hf+hb). (2.28)

    Equations (2.25) and (2.26) can be rewritten as

    2ux2|i=δ2xuihfhb33ux3|i112h3f+h3bhf+hb4ux4|i160(h2f+h2b)(hf+hb)5ux5|i+O(h5f+h5bhf+hb), (2.29)

    and

    ux|i=δxuihfhb63ux3|ihfhb(hfhb)244ux4|i+O((h3f+h3b)hfhbhf+hb). (2.30)

    Substituting Eqs (2.29) and (2.30) into the governing equation (1.1), we can get the difference scheme for the elliptic interface equation

    βδ2xuiβxδxui+kiui+τ=fi, (2.31)

    where

    τ=Q13ux3|i+Q24ux4|i+Q35ux5|i+O(h4), (2.32)

    and

    Q1=βhfhb3βxhfhb6,Q2=β112h3f+h3bhf+hbβxhfhb(hfhb)24,Q3=β160(h2f+h2b)(hf+hb).

    To obtain the higher-order accuracy discrete scheme, it is necessary to handle the third- and fourth-order derivatives of the term τ. Substituting (2.4) and (2.5) into (2.31) and rearranging it, we have

    Aiδ2xui+Bδxui+Ciui=Fi+O((h2f+h2b)(hf+hb)), (2.33)

    where

    Ai=β+Q1ϕ1+Q2ψ1,Bi=βx+Q1ϕ2+Q2ψ2,Ci=ci+Q1ϕ3+Q2ψ3,Fi=fiQ1(ϕ4fi+ϕ5fx,i)Q2ψ4.

    For the irregular points k and k+1, we use the interface α as one of the points in the compact stencil as shown in Figure 5. Thus, the 3-point compact stencil for nodes k and k+1 are k1,k,α and α+,k+1,k+2, respectively. From the above analysis, the difference operators for the first- and second-order derivatives for nodes k and k+1 are respectively as follows

    δxuk=h2buh2fuk1+(h2fh2b)ukh2bhf+h2fhb,δ2xuk=2hbu+hfuk1(hf+hb)ukhfhb(hf+hb), (2.34)
    Figure 5.  Stencil of interface unfitted with a body.

    and

    δxuk+1=h2buk+2h2fu++(h2fh2b)uk+1h2bhf+h2fhb,δ2xuk+1=2hbuk+2+hfu+(hf+hb)uk+1hfhb(hf+hb). (2.35)

    Substituting (2.34) and (2.35) into (2.33), we can get the higher-order compact scheme for the irregular point k:

    ˉAkuk1+ˉBkuk+ˉCku=ˉFk+O((h2f+h2b)(hf+hb)), (2.36)

    where

    ˉAk=2hb(hf+hb)Akh2fh2bhf+h2fhbBk,ˉBk=2hfhbAk+h2fh2bh2bhf+h2fhbBk+Ck,ˉCk=2hf(hf+hb)Ak+h2bh2bhf+h2fhbBk,ˉFk=Fk.

    Similarly, the higher-order compact scheme for the irregular point k+1 is as follows:

    ˉAk+1u++ˉBk+1uk+1+ˉCk+1uk+2=ˉFk+1+O((h2f+h2b)(hf+hb)), (2.37)

    where

    ˉAk+1=2hb(hf+hb)Ak+1h2fh2bhf+h2fhbBk+1,ˉBk+1=2hfhbAk+1+h2fh2bh2bhf+h2fhbBk+1+Ck+1,ˉCk+1=2hf(hf+hb)Ak+1+h2bh2bhf+h2fhbBk+1,ˉFk+1=Fk+1.

    From the above analysis, we have constructed the higher-order compact scheme at regular points and irregular points. It should be noted that there are two auxiliary unknown quantities u+ and u in the scheme for the irregular points k and k+1. Therefore, we have to find a way to deal with the interface connection conditions, so as to obtain a high-precision approximation for the auxiliary unknowns u+ and u on both sides of the interface.

    The first order derivative on both sides of the interface, ux|+ and ux|, can respectively be given as

    ux|+=uk+1u+hr1hr12!2ux2|+h2r13!3ux3|+h3r14!4ux4|+h4r15!5ux5|++O(h5r), (3.1)

    and

    ux|=uukhl1+hl12!2ux2|h2l13!3ux3|+h3l14!4ux4|h4l15!5ux5|+O(h5l), (3.2)

    where hl1=αxk, hr1=xk+1α, hl2=hl1+h and hr2=hr1+h, as shown in Figure 6.

    Figure 6.  Interface and the irregular nodes.

    Further, the second-order derivative on both sides of the interface, 2ux2|+ and 2ux2|, can respectively be given as

    2ux2|+=2hr1hr2h(hr2uk+1hr1uk+2hu+)hr1+hr233ux3|+h2r1+hr1hr2+h2r2124ux4|++O(h3r1+h2r2hr1+h2r1hr2+h3r2), (3.3)

    and

    2ux2|=2hl1hl2h(hl2ukhl1uk1hu)+hl1+hl233ux3|h2l1+hl1hl2+h2l2124ux4|+O(h3l1+h2l2hl1+h2l1hl2+h3l2). (3.4)

    Substituting (3.3) and (3.4) into (3.1) and (3.2) respectively, we can get

    ux|+={(1hr1+1hr2)u++(1hr1+1h)uk+1hr1hr2huk+2}+hr1hr23!3ux3|++h2r1hr2+h2r2hr14!4ux4|++O(h3r1+h2r2hr1+h2r1hr2+h3r2), (3.5)

    and

    ux|={(1hl1+1hl2)u(1hl1+1h)uk+hl1hl2huk1}+hl1hl23!3ux3|+h2l1hl2+h2l2hl14!4ux4|++O(h3l1+h2l2hl1+h2l1hl2+h3l2). (3.6)

    Keep the first term in (3.6) and (3.5) and discard the other terms. We can get the approximation of the first order derivative on both sides of the interface with second-order accuracy:

    ux=ρl,1u+ρl,2uk+ρl,3uk1+O(h2), (3.7)

    and

    u+x=ρr,1u++ρr,2uk+1+ρr,3uk+2+O(h2), (3.8)

    where

    ρl,1=1hl1+1hl2,ρl,2=(1hl1+1h),ρl,3=hl1hhl2,ρr,1=(1hr1+1hr2),ρr,2=1hr1+1h,ρr,3=hr1hhr2.

    Keeping the first two terms in (3.6) and (3.5) and discarding the other terms, and according to the (2.16), we can get the approximation of the first derivative with third order accuracy

    ux=ˉρl,1u+ˉρl,2uk+ˉρl,3uk1+ˉρl,4+O(h3), (3.9)

    and

    u+x=ˉρr,1u++ˉρr,2uk+1+ˉρr,3uk+2+ˉρr,4+O(h3), (3.10)

    where

    ˉρl,1=(ρ1,1+hl1hl26Φ2)/ˉDl,ˉρl,2=ρ1,2/ˉDl,ˉρl,3=ρ1,3/ˉDl,ˉρl,4=hl1hl26Φ3/ˉDl,ˉDl=1hl1hl26Φ1,ˉρr,1=(ρr,1+hr1hr26Φ+2)/ˉDr,ˉρr,2=ρr,2/ˉDr,ˉρr,3=ρr,3/ˉDr,ˉρl,4=hr1hr26Φ+3/ˉDr,ˉDr=1hr1hr26Φ+1.

    Following the same manner as in the above subsection, we can get the approximation of the first derivative with fourth-order accuracy

    ux=ˉˉρl,1u+ˉˉρl,2uk+ˉˉρl,3uk1+ˉˉρl,4+O(h4), (3.11)

    and

    u+x=ˉˉρr,1u++ˉˉρr,2uk+1+ˉˉρr,3uk+2+ˉˉρr,4+O(h4), (3.12)

    where

    ˉˉρl,1=(ρl,1+hl1hl26Φ2+h2l2hl1+h2l1hl224Ψ2)/ˉˉDl,ˉˉρl,2=ρl,2/ˉˉDl,ˉˉρl,3=ρl,3/ˉˉDl,ˉˉρl,4=(hl1hl26Φ3+h2l2hl1+h2l1hl224Ψ3)/ˉˉDl,ˉˉρr,1=(ρr,1hl1hl26Φ1h2l2hl1+h2l1hl224Ψ1)/ˉˉDr,ˉˉρr,2=ρr,2/ˉˉDr,ˉˉρr,3=ρr,3/ˉˉDr,ˉˉρr,4=(hr1hr26Φ+3h2r2hr1+h2r1hr224Ψ+3)/ˉˉDr,ˉˉDl=1hl1hl26Φ1h2r2hr1+h2r1hr224Ψ+1,ˉˉDr=1hr1hr26Φ+1h2r2hr1+h2r1hr224Ψ+1.

    The implicit connecting condition on the interface is given as

    {[u]=λβux,β+u+x=βux. (4.1)

    For the interface fitted mesh, the approximation of the first-order derivative on both sides of the interface can respectively be written as

    ux|=ρl,1u+ρl,2uk1+ρl,3+O(hml), (4.2)
    ux|+=ρr,1u++ρr,2uk+1+ρr,3+O(hmr). (4.3)

    Discarding the terms O(hml) and O(hmr) in (4.2) and (4.3), and substituting them into (4.1), the discrete format of the implicit connecting condition on the interface can be given as

    {u+u=λβ(ρl,1u+ρl,2uk1+ρl,3),β+(ρr,1u++ρr,2uk+1+ρr,3)=β(ρl,1u+ρl,2uk1+ρl,3).

    For brevity, rearrange the above equation and establish linear equations with u+ and u as variables as follows:

    {a11u+a12u=b1,a21u++a22u=b2,

    where

    a11=1,a12=1+λβρl,1,a21=β+ρr,1,a22=βρl,1.b1=λβ(ρl,2uk1+ρl,3),b2=β(ρl,2uk1+ρl,3)β+(ρr,2uk+1+ρr,3).

    Solving the above linear system, we can get the approximate of each of the unknowns u+ and u on both sides of the interface with m-th order accuracy:

    u+=τ+1uk1+τ+2uk+1+τ+3, (4.4)
    u=τ1uk1+τ2uk+1+τ3, (4.5)

    where

    τ+1=(a22λ+a12)βρl,2a11a22+a12a21,τ+2=a12β+ρr,2a11a22+a12a21,τ+3=(a22λ+a12)βρl,3a12β+ρr,3a11a22+a12a21,τ1=(a11a21λ)βρl,2a11a22+a12a21,τ2=a11β+ρr,2a11a22+a12a21,τ3=(a11a21λ)βρl,3a11β+ρr,3a11a22+a12a21.

    For the interface-cut mesh, the approximation of the first-order derivative on both sides of the interface can be written as

    ux|=ρl,1u+ρl,2uk+ρl,3uk1+ρl,4+O(hml), (4.6)
    ux|+=ρr,1u++ρr,2uk+1+ρr,3uk+2+ρr,4+O(hmr). (4.7)

    In a similar manner as above, the approximate of each of the unknowns u+ and u on both sides of the interface with m-th-order accuracy is given as

    u+=τ+1uk1+τ+2uk+τ+3uk+1+τ+4uk+2+τ+5, (4.8)
    u=τ1uk1+τ2uk+τ3uk+1+τ4uk+2+τ5, (4.9)

    where

    τ+1=(a22λ+a12)βρl,3a11a22+a12a21,τ+2=(a22λ+a12)βρl,2a11a22+a12a21,τ+3=a12β+ρr,2a11a22+a12a21,τ+4=a12β+ρr,3a11a22+a12a21,τ+5=a22λβρl,4+a12(βρl,4β+ρr,4)a11a22+a12a21,τ1=(a11a21λ)βρl,3a11a22+a12a21,τ2=(a11a21λ)βρl,2a11a22+a12a21,τ3=a11β+ρr,2a11a22+a12a21,τ4=a11β+ρr,3a11a22+a12a21,τ5=a11(βρl,4β+ρr,4)a21λβρl,4a11a22+a12a21.

    and

    a11=1,a12=1+λβρl,1,a21=β+ρr,1,a22=βρl,1.b1=λβ(ρl,2uk+ρl,3uk1+ρl,4),b2=β(ρl,2uk+ρl,3uk1+ρl,4)β+(ρr,2uk+1+ρr,3uk+2+ρr,4).

    From the above sections, we establish the fourth-order compact finite-difference scheme for non-irregular mesh nodes:

    Aiui1+Biui+Ciui+1=Fi, (5.1)

    for

    i{{k1,k+1},forbodyfittedmesh,{k,k+1},fornonbodyfittedmesh.

    In the interface-fitted mesh case, the higher-order compact schemes at the irregular mesh nodes k1 and k+1 are given as (2.7) and (2.8), respectively. Eliminating the auxiliary qualities u+ and u by substituting the formulas (4.4) and (4.5) into (2.7) and (2.8). The schemes for the irregular mesh nodes are given as

    Ak1uk2+Bk1uk1+Ck1uk+1=Fk1, (5.2)

    and

    Ak+1uk1+Bk+1uk+1+Ck+1uk+2=Fk+1, (5.3)

    where

    Ak1=ˉAk1,Bk1=ˉBk1+ˉCk1τ1,Ck1=ˉCk1τ2,Ak+1=ˉAk+1τ+1,Bk+1=ˉBk+1+ˉAk+1τ+2,Ck+1=ˉCk+1,Fk1=ˉFk1ˉCk1τ3,Fk+1=ˉFk+1ˉAk+1τ+3.

    Similarly, in the non-body-fitted case, the higher-order compact schemes at the irregular mesh nodes k and k+1 are respectively given as (2.36) and (2.37). Substitute the formulas (4.8) and (4.9) into (2.36) and (2.37) to eliminate the auxiliary qualities u+ and u. The higher-order compact schemes at the irregular mesh nodes k and k+1 are respectively given as

    Akuk1+Bkuk+Ckuk+1+Dkuk+2=Fk, (5.4)

    and

    Ak+1uk1+Bk+1uk+Ck+1uk+1+Dk+1uk+2=Fk+1, (5.5)

    where

    Ak=ˉAk+ˉCτ1,Bk=ˉBk+ˉCτ2,Ck=ˉCkτ3,Dk=ˉCkτ4,Ak+1=ˉAk+1τ+1,Bk+1=ˉAk+1τ+2,Ck+1=ˉBk+1+ˉAk+1τ+3,Dk+1=ˉCk+1+ˉAk+1τ+4,Fk=ˉFkˉCτ5,Fk+1=ˉFk+1ˉAk+1τ+5.

    Obviously, the resulting discrete linear equations have a tri-diagonal form for the body-fitted case and block-diagonal form for the non-body-fitted case, which can be efficiently solved by using existing numerical methods, such as the forward and backward sweep method, BiCGstab method, etc.

    In this section, we use several numerical experiments to demonstrate the performance of the discrete schemes.

    Example 1. Consider the computational domain Ω=[0,1], and the solution is separated into two parts by the interface at x=α, where α=0.5 and α=0.52323, for the interface-fitted mesh and interface-cut mesh cases, respectively. The analytical solution of this problem is given by

    u(x,y)={ex2,x(0,α),κex2,x(α,1).

    The diffusion coefficient is defined as follows

    β={κ,x(0,α),1,x(α,1).

    Conservation of the flux on the interface satisfies

    [βux]=β+u+xβux=2κxex22κxex2=0,atx=α.

    The coefficient λ is given as

    λ=κeα2α22καeα2.

    Tables 1 and 2 compare the L2 and L errors in second-, third- and fourth-order formats on fitted and non-fitted mesh, respectively. The tables show that the numerical results of the three methods achieve theoretical accuracy on both fitted and non-fitted mesh.

    Table 1.  Comparison of L2 and L errors in second-, third- and fourth-order formats on fitted mesh for Example 1, κ = 100.
    Mesh Second-order scheme Third-order scheme Fourth-order scheme
    L2 Rate L Rate L2 Rate L Rate L2 Rate L Rate
    11 4.8e-1 8.7e-1 2.2e-2 4.7e-2 2.4e-3 4.5e-3
    21 1.2e-1 2.2 2.2e-1 2.1 2.4e-3 3.2 5.6e-3 3.1 1.4e-4 4.1 2.8e-4 4.0
    41 2.8e-2 2.1 5.5e-2 2.0 2.8e-4 3.1 6.8e-4 3.0 8.8e-6 4.0 1.7e-5 4.0
    81 6.9e-3 2.0 1.4e-2 2.0 3.4e-5 3.1 8.4e-5 3.0 5.4e-7 4.0 1.1e-6 4.0
    161 1.7e-3 2.0 3.4e-3 2.0 4.2e-6 3.0 1.0e-5 3.0 3.3e-8 4.0 6.8e-8 4.0

     | Show Table
    DownLoad: CSV
    Table 2.  Comparison of L2 and L errors in second-, third- and fourth-order formats on non-fitted mesh Example 1, κ = 100.
    Mesh Second-order scheme Third-order scheme Fourth-order scheme
    L2 Rate L Rate L2 Rate L Rate L2 Rate L Rate
    11 8.9e-2 1.8e-1 9.8e-2 2.7e-1 2.4e-3 4.5e-3
    21 6.9e-2 0.4 1.2e-1 0.6 2.4e-2 2.0 4.5e-2 2.6 1.5e-4 4.0 2.8e-4 4.0
    41 2.8e-2 1.3 5.5e-2 1.1 3.5e-3 2.8 7.3e-3 2.6 8.8e-6 4.0 1.8e-5 4.0
    81 6.8e-3 2.2 1.3e-2 2.0 4.3e-4 3.1 8.4e-4 3.1 5.4e-7 4.0 1.1e-6 4.0
    161 1.5e-3 2.1 3.0e-3 2.2 5.5e-5 2.9 1.0e-4 3.0 3.4e-8 4.0 6.8e-8 5.0

     | Show Table
    DownLoad: CSV

    Due to the strong discontinuity of physical quantities at the interface, it is crucial for the numerical scheme to be stable and robust. Tables 3 and 4 compare the errors of u+ and u for different methods on fitted and non-fitted mesh, respectively. It can be seen from the tables that the physical quantities on both sides of the interface differ in magnitude by a factor of 100, but the numerical scheme proposed can capture the discontinuity of physical quantities on both sides of the interface with a high degree of accuracy.

    Table 3.  Comparison of the errors of u+ and u for different methods on fitted mesh for Example 1, κ = 100.
    (u,u+) = (1.284025,128.4025)
    Mesh 11 21 41 81 161
    Second e (8.6e-4, 8.7e-1) (2.1e-4, 2.2e-1) (5.3e-5, 5.5e-2) (1.3e-5, 1.4e-2) (3.3e-6, 3.4e-3)
    e2h/eh - (4.0, 4.0) (4.0, 4.0) (4.0, 4.0) (4.0, 4.0)
    Third e (8.6e-4, 4.7e-2) (1.1e-4, 5.6e-3) (1.4e-5, 6.8e-4) (1.7e-6, 8.4e-5) (2.1e-7, 1.0e-5)
    e2h/eh - (8.0, 8.4) (7.9, 8.2) (8.0, 8.1) (8.0, 8.1)
    Fourth e (1.2e-5, 4.5e-3) (7.3e-7, 2.8e-4) (4.6e-8, 1.8e-5) (2.9e-9, 1.1e-6) (1.8e-10, 6.8e-8)
    e2h/eh - (15.9, 16.0) (16.0, 16.0) (16.0, 16.1) (16.1, 16.0)

     | Show Table
    DownLoad: CSV
    Table 4.  Comparison of the errors of u+ and u for different methods on non-fitted mesh for Example 1, κ = 100.
    (u,u+) = (1.284025,128.4025)
    mesh 11 21 41 81 161
    Second e (4.7e-3, 1.3e-1) (1.6e-3, 1.1e-1) (1.3e-3, 5.5e-2) (2.7e-4, 1.3e-2) (6.4e-5, 3.0e-3)
    e2h/eh - (3.6, 4.2) (3.8, 4.1) (3.9, 4.1) (4.2, 4.5)
    Third e (8.7e-4, 4.7e-2) (1.1e-4, 5.6e-3) (1.4e-5, 6.8e-4) (1.7e-6, 8.4e-5) (2.1e-7, 1.0e-5)
    e2h/eh - (8.0, 8.4) (7.9, 8.2) (8.0, 8.1) (8.0, 8.1)
    Fourth e (1.2e-5, 4.5e-3) (7.3e-7, 2.8e-4) (4.6e-8, 1.8e-5) (2.9e-9, 1.1e-6) (1.8e-10, 6.8e-8)
    e2h/eh - (15.9, 16.0) (16.0, 16.0) (16.0, 16.1) (16.1, 16.0)

     | Show Table
    DownLoad: CSV

    Figure 7 compares the exact and numerical solutions on fitted and non-fitted mesh for Problem 1 at k=100 with a number of intervals of 32. It is shown that the numerical solutions are well matched with the exact solutions, in spite of the fact that the amount of physical jumps on both sides of the interface are much larger. Figures 8 and 9 errors by the fourth-order, second-order and third-order scheme on fitted and non-fitted mesh, respectively. It can be seen that the errors are decrease as the number of mesh points increases.

    Figure 7.  Comparison of exact and numerical solutions on fitted and non-fitted mesh for Problem 1.
    Figure 8.  Comparison of errors by the fourth-order scheme on fitted mesh and non-fitted mesh for Example 1.
    Figure 9.  Comparison of errors by the second-order and third-order schemes for Example 1.

    Example 2. Consider the elliptic interface problem with a variable diffusion coefficient, and the solution is separated into two parts by the interface at x=α, where α=0.5 and α=0.53232, for the interface-fitted mesh and interface-cut mesh cases, respectively. The analytical solution of this problem is given by

    u(x,y)={x2ex,x(0,α),κex,x(α,1).

    The diffusion coefficient is defined as follows

    β={κ,x(0,α),2x+x2,x(α,1).

    Conservation of the flux on the interface satisfies

    [βux]=β+u+xβux=κ(2x+x2)exκ(2x+x2)ex=0,atx=α.

    The coefficient λ is given as

    λ=κα2κ(2α+α2).

    Tables 5 and 6 compare the L2 and L errors in second-, third- and fourth-order formats on fitted and non-fitted mesh, respectively. It can be seen that the numerical format is able to achieve theoretical accuracy on fitted mesh. The numerical accuracy on the non-fitted mesh is slightly lower than the theoretical accuracy, probably due to the unbalanced grid step on both sides of the interface, which leads to the instability of the format.

    Table 5.  Comparison of L2 and L errors in second-, third- and fourth-order formats on fitted mesh for Example 2, κ = 100.
    Mesh Second-order scheme Third-order scheme Fourth-order scheme
    L2 Rate L Rate L2 Rate L Rate L2 Rate L Rate
    21 1.4e-2 2.9e-2 9.9e-5 2.4e-4 3.2e-6 6.1e-6
    41 3.4e-3 2.0 7.3e-3 2.0 1.1e-5 3.1 2.9e-5 3.1 2.0e-7 4.0 4.0e-7 4.0
    81 8.3e-4 2.0 1.8e-3 2.0 1.4e-6 3.1 3.5e-6 3.0 1.2e-8 4.0 2.5e-8 4.0
    161 2.1e-4 2.0 4.6e-4 2.0 1.7e-7 3.0 4.3e-7 3.0 7.5e-10 4.0 1.5e-9 4.1
    321 5.1e-5 2.0 1.1e-4 2.0 2.1e-8 3.0 5.4e-8 3.0 4.6e-11 4.0 9.2e-11 4.1

     | Show Table
    DownLoad: CSV
    Table 6.  Comparison of L2 and L errors in second-, third- and fourth-order formats on non-fitted mesh for Example 2, κ = 100.
    Mesh Second-order scheme Third-order scheme Fourth-order scheme
    L2 Rate L Rate L2 Rate L Rate L2 Rate L Rate
    21 8.9e-3 1.5e-2 1.3e-3 5.9e-3 1.1e-3 5.2e3
    41 3.6e-3 1.3 7.6e-3 1.0 3.2e-5 5.3 7.5e-5 6.3 3.2e-5 5.1 8.8e-5 5.9
    81 8.3e-4 2.1 1.8e-3 2.1 3.7e-6 3.1 1.3e-5 2.5 3.8e-6 3.1 1.5e-5 2.6
    161 1.8e-4 2.2 4.0e-4 2.2 4.2e-7 3.1 3.8e-6 1.8 4.8e-7 3.0 4.0e-6 1.9
    321 2.8e-5 2.7 5.2e-5 2.9 9.5e-8 2.2 1.7e-6 1.2 9.9e-8 2.3 1.7e-6 1.3

     | Show Table
    DownLoad: CSV

    Figure 10 compares the exact and numerical solutions on fitted and non-fitted mesh for Problem 1 at k=100 with a number of intervals of 41. It is shown that the numerical solutions are well matched with the exact solutions, in spite of the fact that the amount of physical jumps on both sides of the interface are much larger. Figure 11 compares errors by the fourth-order and third-order schemes on fitted and non-fitted mesh, respectively. As seen in the table, the errors in the fourth-order scheme are significantly lower than those in the third-order scheme.

    Figure 10.  Comparison of exact and numerical solutions for Problem 2.
    Figure 11.  Comparison of errors under different kinds schemes for Problem 2.

    Table 7 compares the length ratio of the interface-cut mesh interval. The table shows that as the number of grids changes, the spacing ratio changes significantly and affects the accuracy of the numerical scheme.

    Table 7.  Comparison of the length ratio hf/hb of interface-cut mesh.
    Mesh 21 41 81 161 321
    hf/hb 1.48 2.80 1.72e-2 1.36 2.30

     | Show Table
    DownLoad: CSV

    Table 8 compares the condition numbers of the linear system matrix A for the fitted and non-fitted meshes. It can be seen that the condition number of the matrix does not differ much for the same grid number, so the proposed numerical scheme is stable.

    Table 8.  Comparison of condition numbers of the linear system on fitted and non-fitted mesh.
    Mesh 21 41 81 161 321
    Fitted 2098 6746 20081 63749 147695
    Non-fitted 4108 14355 31339 107609 341350

     | Show Table
    DownLoad: CSV

    Example 3. Consider the elliptic interface problem with a variable diffusion coefficient, and the solution is separated into two parts by the interface at x=α, where α=0.5 and α=0.53232, for the interface-fitted mesh and interface-cut mesh cases, respectively. The analytical solution of this problem is given by

    u(x,y)={sinxex2,x(0,α),κex2,x(α,1).

    The diffusion coefficient is defined as follows

    β={2κx,x(0,α),cosx+2xsinx,x(α,1).

    Conservation of the flux on the interface satisfies

    [βux]=β+u+xβux=2κx(2xsinx+cosx)ex22κx(2xsinx+cosx)ex2=0.

    The coefficient λ is given as

    λ=κsinα2κα(cosα+2αsinα).

    Tables 9 and 10 compare the L2 and L errors in second-, third- and fourth-order formats on fitted and non-fitted mesh, respectively. It is shows that the numerical format is able to achieve theoretical accuracy on fitted mesh.

    Table 9.  Comparison of L2 and Linfty errors in second-, third- and fourth-order formats on fitted mesh for Example 3, κ = 1000.
    Mesh Second-order scheme Third-order scheme Fourth-order scheme
    L2 Rate L Rate L2 Rate L Rate L2 Rate L Rate
    21 1.4 2.6 2.1e-2 5.0e-2 1.8e-3 3.4e-3
    41 3.4e-1 2.0 6.4e-1 2.0 2.4e-3 3.1 6.0e-3 3.1 1.1e-4 4.0 2.1e-4 4.0
    81 8.4e-2 2.0 1.6e-1 2.0 2.9e-4 3.1 7.4e-4 3.0 6.9e-6 4.0 1.3e-5 4.0
    161 2.1e-2 2.0 4.0e-2 2.0 3.5e-5 3.0 9.1e-5 3.0 4.4e-7 4.0 8.5e-7 4.0
    321 5.2e-3 2.0 1.0e-2 2.0 4.3e-6 3.0 1.1e-5 3.0 2.9e-8 3.9 5.5e-8 3.9

     | Show Table
    DownLoad: CSV
    Table 10.  Comparison of L2 and Linfty errors in second-, third- and fourth-order formats on non-fitted mesh for Example 3, κ = 1000.
    Mesh Second-order scheme Third-order scheme Fourth-order scheme
    L2 Rate L Rate L2 Rate L Rate L2 Rate L Rate
    21 9.9e-1 1.7 4.0e-2 1.4e-1 1.4e-2 6.3e-2
    41 3.4e-1 1.6 6.4e-1 1.4 8.2e-4 5.6 1.9e-3 6.2 1.1e-3 3.7 2.5e-3 4.6
    81 8.1e-2 2.1 1.6e-1 2.0 5.3e-5 4.0 1.2e-4 4.0 1.1e-4 3.2 2.9e-4 3.1
    161 1.9e-2 2.1 3.6e-2 2.1 1.0e-5 2.4 5.9e-5 1.0 1.1e-5 3.3 3.0e-5 3.3
    321 3.6e-3 2.4 6.3e-3 2.5 1.8e-6 2.5 3.9e-5 0.6 1.2e-6 3.3 1.7e-6 4.1

     | Show Table
    DownLoad: CSV

    Figure 12 compares the exact and numerical solutions on fitted and non-fitted mesh for Problem 1 at k=100 with a number of intervals of 41. It can be seen that the error in matching the numerical solution to the exact solution is essentially small despite the large physical jumps on both sides of the interface. Figure 13 compares errors by the fourth-order and third-order schemes on fitted and non-fitted mesh, respectively. As can be seen in the table, the numerical accuracy of the fourth-order scheme is better than that of the third-order scheme.

    Figure 12.  Comparison of exact and numerical solutions for Problem 3.
    Figure 13.  Errors under the third-order scheme (a) and fourth-order scheme (b) for Problem 3, with different mesh sizes.

    The 1D elliptic interface problem with imperfect contact is characterized by the fact that the jump quantity of the solution is unknown and is related to the flux across the interface. In this paper, a class of higher-order finite-difference schemes is constructed for interface fitted and non-fitted meshes, respectively. The second-, third-, and fourth-order approximations of the jump conditions are provided by using the connected jump conditions and their higher-order derivatives. Some numerical experiments were carried out to illustrate the accuracy and stability of the present method. The numerical results show that the proposed scheme is able to capture the solution of the discontinuous case of the interface body with theoretical accuracy. The theoretical accuracy has been achieved for the elliptic interface problem with implicit interface connection conditions.

    This work was partly supported by the National Natural Science Foundation (12161067, 12261067, 12001015, 62201298, 51961031), the National Natural Science Foundation of China Youth Fund Project (11801287), the Inner Mongolia Autonomous Region Youth Science and Technology Talents Support Program (NJYT20B15), the Inner Mongolia Scientific Fund Project (2020MS06010, 2021LHMS01006, 2022MS01008), and Innovation Fund Project of Inner Mongolia University of Science and Technology-Excellent Youth Science Fund Project (2019YQL02).

    The authors declare that the publication of this paper does coincide with any confilict of interest.



    [1] A. Ouardghi, M. El-Amrani, M. Seaid, An enriched Galerkin-characteristics finite element method for convection-dominated and transport problems, Appl. Numer. Math., 167 (2021), 119–142. https://doi.org/10.1016/j.apnum.2021.04.018 doi: 10.1016/j.apnum.2021.04.018
    [2] D. Broersen, R. Stevenson, A robust Petrov-Galerkin discretization of convection–diffusion equations, Comput. Math. Appl., 68 (2014), 1605–1618. https://doi.org/10.1016/j.camwa.2014.06.019 doi: 10.1016/j.camwa.2014.06.019
    [3] A. Cangiani, E. H. Georgoulis, S. Giani, S. Metcalfe, hp-adaptive discontinuous Galerkin methods for non-stationary convection–diffusion problems, Comput. Math. Appl., 78 (2019), 3090–3104. https://doi.org/10.1016/j.camwa.2019.04.002 doi: 10.1016/j.camwa.2019.04.002
    [4] A. El Kacimi, O. Laghrouche, Numerical modelling of elastic wave scattering in frequency domain by the partition of unity finite element method, Int. J. Numer. Methods Eng., 77 (2009), 1646–1669. https://doi.org/10.1002/nme.2471 doi: 10.1002/nme.2471
    [5] X. Xiao, X. Feng, Z. Li, A gradient recovery-based adaptive finite element method for convection-diffusion-reaction equations on surfaces, Int. J. Numer. Methods Eng., 120 (2019), 901–917. https://doi.org/10.1002/nme.6163 doi: 10.1002/nme.6163
    [6] K. Oldham, J. Spanier, The Fractional Calculus Theory and Applications of Differentiation and Integration to Arbitrary Order, Amsterdam: Elsevier, 1974.
    [7] M. Caputo, M. Fabrizio, A new definition of fractional derivative without singular kernel, Progr. Fract. Differ. Appl., 1 (2015), 73–85. http://dx.doi.org/10.12785/pfda/010201 doi: 10.12785/pfda/010201
    [8] A. Atangana, D. Baleanu, New fractional derivatives with non-local and non-singular kernel: theory and applications to heat transfer model, Therm. Sci., 20 (2016), 763–769.
    [9] J. Hadamard, Essai sur l'etude des fonctions donnes par leur developpment de Taylor, J. Pure Appl. Math., 4 (1892), 101–186.
    [10] J. Sabatier, O. P. Agrawal, J. A. T. Machado, Advances in Fractional Calculus, Dordrecht: Springer, 2007.
    [11] A. D. Freed, K. Diethelm, Fractional calculus in biomechanics: a 3D viscoelastic model using regularized fractional derivative kernels with application to the human calcaneal fat pad, Biomech. Model. Mechanobiol., 5 (2006), 203–215. https://doi.org/10.1007/s10237-005-0011-0 doi: 10.1007/s10237-005-0011-0
    [12] M. M. Meerschaert, E. Scalas, Coupled continuous time random walks in finance, Phys. A Stat. Mech. Appl., 370 (2006), 114–118. https://doi.org/10.1016/j.physa.2006.04.034 doi: 10.1016/j.physa.2006.04.034
    [13] O. Sadek, L. Sadek, S. Touhtouh, A. Hajjaji, The mathematical fractional modeling of TiO-2 nanopowder synthesis by sol-gel method at low temperature, Math. Model. Comput., 9 (2022), 616–626. https://doi.org/10.23939/mmc2022.03.616 doi: 10.23939/mmc2022.03.616
    [14] W. M. Ahmad, R. El-Khazali, Fractional-order dynamical models of love, Chaos Solitons Fract., 33 (2007), 1367–1375. https://doi.org/10.1016/j.chaos.2006.01.098 doi: 10.1016/j.chaos.2006.01.098
    [15] A. Atangana, D. Baleanu, New fractional derivatives with nonlocal and non-singular kernel: theory and application to heat transfer model, Therm. Sci., 20 (2016), 757–763.
    [16] F. Gao, X. J. Yang, Fractional Maxwell fluid with fractional derivative without singular kernel, Therm. Sci., 20 (2016), 871–877.
    [17] J. Losada, J. J. Nieto, Properties of a new fractional derivative without singular kernel, Progr. Fract. Differ. Appl., 1 (2015), 87–92. http://dx.doi.org/10.12785/pfda/010202 doi: 10.12785/pfda/010202
    [18] X. J. Yang, F. Gao, J. A. Tenreiro Machado, D. Baleanu, A new fractional derivative involving the normalized sinc function without singular kernel, Eur. Phys. J. Spec. Top., 226 (2017), 3567–3575. https://doi.org/10.1140/epjst/e2018-00020-2 doi: 10.1140/epjst/e2018-00020-2
    [19] T. Abdeljawad, D. Baleanu, Integration by parts and its applications of a new nonlocal fractional derivative with Mittag-Leffler nonsingular kernel, J. Nonlinear Sci. Appl., 10 (2017), 1098–1107.
    [20] T. Abdeljawad, D. Baleanu, Monotonicity results for fractional difference operators with discrete exponential kernels, Adv. Differ. Equ., 2017 (2017), 78. https://doi.org/10.1186/s13662-017-1126-1 doi: 10.1186/s13662-017-1126-1
    [21] T. Abdeljawad, D. Baleanu, On fractional derivatives with exponential kernel and their discrete versions, Rep. Math. Phys., 80 (2017), 11–27. https://doi.org/10.1016/S0034-4877(17)30059-9 doi: 10.1016/S0034-4877(17)30059-9
    [22] L. Sadek, A cotangent fractional derivative with the application, Fractal Fract., 7 (2023), 444. https://doi.org/10.3390/fractalfract7060444 doi: 10.3390/fractalfract7060444
    [23] L. Sadek, A. S. Bataineh, H. Talibi Alaoui, I. Hashim, The novel Mittag-Leffler–Galerkin method: application to a riccati differential equation of fractional order, Fractal Fract., 7 (2023), 302. https://doi.org/10.3390/fractalfract7040302 doi: 10.3390/fractalfract7040302
    [24] R. Khalil, M. Al Horani, A. Yousef, M. Sababheh, A new definition of fractional derivative, J. Comput. Appl. Math., 264 (2014), 65–70. https://doi.org/10.1016/j.cam.2014.01.002 doi: 10.1016/j.cam.2014.01.002
    [25] T. Abdeljawad, On conformable fractional calculus, J. Comput. Appl. Math., 279 (2015), 57–66. https://doi.org/10.1016/j.cam.2014.10.016 doi: 10.1016/j.cam.2014.10.016
    [26] A. Atangana, D. Baleanu, A. Alsaedi, New properties of conformable derivative, Open Math., 13 (2015), 889–898. https://doi.org/10.1515/math-2015-0081 doi: 10.1515/math-2015-0081
    [27] O. Naifar, G. Rebiai, A. B. Makhlouf, M. A. Hammami, A. Guezane-Lakoud, Stability analysis of conformable fractional-order nonlinear systems depending on a parameter, J. Appl. Anal., 26 (2020), 287–296. https://doi.org/10.1515/jaa-2020-2025 doi: 10.1515/jaa-2020-2025
    [28] A. Kütahyalioglu, F. Karakoç, Exponential stability of Hopfield neural networks with conformable fractional derivative, Neurocomputing, 456 (2021), 263–267. https://doi.org/10.1016/j.neucom.2021.05.076 doi: 10.1016/j.neucom.2021.05.076
    [29] Z. Hammouch, R. R. Rasul, A. Ouakka, A. Elazzouzi, Mathematical analysis and numerical simulation of the Ebola epidemic disease in the sense of conformable derivative, Chaos Solitons Fract., 158 (2022), 112006. https://doi.org/10.1016/j.chaos.2022.112006 doi: 10.1016/j.chaos.2022.112006
    [30] H. Zhao, T. Li, P. Cui, On stability for conformable fractional linear system, In: 2020 39th Chinese Control Conference, 2020,899–903. https://doi.org/10.23919/CCC50068.2020.9189052
    [31] G. Rebiai, Stability analysis of nonlinear differential equations depending on a parameter with conformable derivative, New Trends Math. Sci., 1 (2021), 44–49. http://dx.doi.org/10.20852/ntmsci.2021.427 doi: 10.20852/ntmsci.2021.427
    [32] L. Sadek, B. Abouzaid, E. M. Sadek, H. T. Alaoui, Controllability, observability and fractional linear-quadratic problem for fractional linear systems with conformable fractional derivatives and some applications, Int. J. Dynam. Control, 11 (2023), 214–228. https://doi.org/10.1007/s40435-022-00977-7 doi: 10.1007/s40435-022-00977-7
    [33] Z. Al-Zhour, Controllability and observability behaviors of a non-homogeneous conformable fractional dynamical system compatible with some electrical applications, Alex. Eng. J., 61 (2022), 1055–1067. https://doi.org/10.1016/j.aej.2021.07.018 doi: 10.1016/j.aej.2021.07.018
    [34] X. Wang, J. Wang, M. Feckan, Controllability of conformable differential systems, Nonlinear Anal. Model. Control, 25 (2020), 658–674. https://doi.org/10.15388/namc.2020.25.18135 doi: 10.15388/namc.2020.25.18135
    [35] J. C. Mayo-Maldonado, G. Fernandez-Anaya, O. F. Ruiz-Martinez, Stability of conformable linear differential systems: a behavioural framework with applications in fractional-order control, IET Control Theory Appl., 14 (2020), 2900–2913. https://doi.org/10.1049/iet-cta.2019.0930 doi: 10.1049/iet-cta.2019.0930
    [36] A. Ben Makhlouf, L. Mchiri, M. Rhaima, M. A. Hammami, Stability of conformable stochastic systems depending on a parameter, Asian J. Control, 25 (2023), 594–603. https://doi.org/10.1002/asjc.2804 doi: 10.1002/asjc.2804
    [37] H. Rezazadeh, H Aminikhah, S. A. Refahi, Stability analysis of conformable fractional systems, Iran. J. Numer. Anal. Optimiz., 7 (2017), 13–32. https://doi.org/10.22067/ijnao.v7i1.46917 doi: 10.22067/ijnao.v7i1.46917
    [38] A. Souahi, A. B. Makhlouf, M. A. Hammami, Stability analysis of conformable fractional-order nonlinear systems, Indagat. Math., 28 (2017), 1265–1274. https://doi.org/10.1016/j.indag.2017.09.009 doi: 10.1016/j.indag.2017.09.009
    [39] Y. Qi, X. Wang, Asymptotical stability analysis of conformable fractional systems, J. Taibah Uni. Sci., 14 (2020), 44–49. https://doi.org/10.1080/16583655.2019.1701390 doi: 10.1080/16583655.2019.1701390
    [40] A. Younus, T. Abdeljawad, T. Gul, On stability criteria of fractal differential systems of conformable type, Fractals, 28 (2020), 2040009. https://doi.org/10.1142/S0218348X20400095 doi: 10.1142/S0218348X20400095
    [41] L. Sadek, Stability of conformable linear infinite-dimensional systems, Int. J. Dynam. Control, 11 (2022), 1276–1284. https://doi.org/10.1007/s40435-022-01061-w doi: 10.1007/s40435-022-01061-w
    [42] M. Yavari, A. Nazemi, On fractional infinite-horizon optimal control problems with a combination of conformable and Caputo–Fabrizio fractional derivatives, ISA Trans., 101 (2020), 78–90. https://doi.org/10.1016/j.isatra.2020.02.011 doi: 10.1016/j.isatra.2020.02.011
    [43] L. Pedram, D. Rostamy, Numerical solutions of the initial boundary value problem for the perturbed conformable time Korteweg-de Vries equation by using the finite element method, Numer. Methods Partial Differ. Equ., 37 (2021), 1449–1463. https://doi.org/10.1002/num.22590 doi: 10.1002/num.22590
    [44] Y. Wang, J. Zhou, Y. Li, Fractional Sobolev's spaces on time scales via conformable fractional calculus and their application to a fractional differential equation on time scales, Adv. Math. Phys., 2016 (2016), 9636491. https://doi.org/10.1155/2016/9636491 doi: 10.1155/2016/9636491
    [45] B. P. Allahverdiev, H. Tuna, Y. Yalçinkaya, Spectral expansion for singular conformable fractional sturm-liouville problem, Math. Commun., 25 (2020), 237–252.
    [46] B. Lucquin, Équations aux dérivées partielles et leurs approximations: niveau M1, Ellipses Éd. Marketing, 2004.
    [47] X. Li, A stabilized element-free Galerkin method for the advection–diffusion–reaction problem, Appl. Math. Lett., 146 (2023), 108831. https://doi.org/10.1016/j.aml.2023.108831 doi: 10.1016/j.aml.2023.108831
    [48] X. Li, Element-free Galerkin analysis of Stokes problems using the reproducing kernel gradient smoothing integration, J. Sci. Comput., 96 (2023), 43. https://doi.org/10.1007/s10915-023-02273-8 doi: 10.1007/s10915-023-02273-8
    [49] X. Li, S. Li, Effect of an efficient numerical integration technique on the element-free Galerkin method, Appl. Numer. Math., 193 (2023), 204–225. https://doi.org/10.1016/j.apnum.2023.07.026 doi: 10.1016/j.apnum.2023.07.026
  • This article has been cited by:

    1. Shougui Zhang, Xiyong Cui, Guihua Xiong, Ruisheng Ran, An Optimal ADMM for Unilateral Obstacle Problems, 2024, 12, 2227-7390, 1901, 10.3390/math12121901
  • Reader Comments
  • © 2023 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(1778) PDF downloads(111) Cited by(4)

Figures and Tables

Figures(4)

Other Articles By Authors

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog