Superconvergence Analysis of a Multiscale Finite Element Method for Elliptic Problems with Rapidly Oscillating Coefficients

Abstract

A new multiscale finite element method is presented for solving the elliptic equations with rapidly oscillating coefficients. The proposed method is based on asymptotic analysis and careful numerical treatments for the boundary corrector terms by virtue of the recovery technique. Under the assumption that the oscillating coefficient is periodic, some superconvergence results are derived, which seem to be never discovered in the previous literature. Finally, some numerical experiments are carried out to demonstrate the efficiency and accuracy of this method, and it is seen that they agree very well with the analytical result.

1. Introduction

In this paper, we consider the following elliptic boundary value problem with rapidly oscillatory coefficients: 𝐿_𝜀𝑢^𝜀≡𝜕𝜕𝑥_𝑖𝑎_𝑖𝑗𝑥𝜀𝜕𝑢^𝜀𝜕𝑥_𝑗(𝑥)𝑢=𝑓,inΩ,_𝜀=𝑔(𝑥),on𝜕Ω, (1) where Ω⊂ℜ² is a smooth-bounded domain, 𝑎_𝑖𝑗(𝜉)∶ℜ²→ℜ is symmetric and satisfies (1)𝜆𝜉²≤𝑎_𝑖𝑗(𝜉)𝜉_𝑖𝜉_𝑗≤𝜆⁻¹𝜉²,∀𝜉∈ℜ²(],(2)𝑎,∃𝜆∈0,1_𝑖𝑗𝜉+𝜉^′=𝑎_𝑖𝑗(𝜉),∀𝜉∈ℜ²,∃𝜉^′∈𝑍²,(3𝑎1≤𝑖,𝑗≤𝑛,)_{𝑖𝑗𝐻1(ℜ2)}≤𝐶,∃𝐶>0, (2) where 𝜉=𝑥/𝜀, 𝜀 is a small scale parameter. This kind of equation has widely been applied in many areas, such as the behavior of flow in porous media or the thermal and mechanical behavior of composite material structure. In practice, the oscillatory coefficients may span many scales to a great extent. In such cases, the direct accurate numerical computation of the solution becomes difficult because it would require a very fine mesh, and it can easily exceed the limit of today’s computer resources because of the requirement of tremendous amount of computer memory and CPU time. Meanwhile, it is desirable to have a numerical method that can solve this equation on a large-scale mesh with capturing the effect of small scales details. Thus, various methods of upscaling or homogenization have been developed.

Based on the homogenization method, there are many discussions [1–4] about the numerical methods of (1). A large amount of examples and applications can also be found in the classical books [5–8], where the formal asymptotic expansions for the limit solution are deduced when 𝜀 is small enough. In these books, the first-order approximation of these expansions is justified by proving sharp error estimates, from which a general method that allowed us to treat some structures with rapidly oscillatory coefficients is also developed. However, the general method cannot effectively compute the boundary corrector on boundary layer. It should be noted that the boundary corrector is the important source of error estimates. In [9], He and Cui present a novel finite element method to solve (1) which can effectively compute the boundary corrector even if the boundary layer is very small. The crucial idea is to combine the numerical approximation of the first-order terms of asymptotic expansions with the numerical approximation of the boundary corrector from different meshes exploiting the need for different levels of resolution. The following result (Theorem  2.13 in [9]) can be obtained.

Lemma 1. Assume that 𝑢^𝜀 is the solution of (1) and 𝑢^{ℎ0,ℎ1,ℎ} is the finite element solution [9]. For all 𝑝,1<𝑝<+∞, there exists a constant 𝐶 such that ∇𝑢^𝜀−𝑢^{ℎ0,ℎ1,ℎ}dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)≤𝐶ℎ₁+ℎ₀+ℎ+𝜀|ln𝜀|^1/2+𝜀^{(2𝑝−1)/2𝑝}×𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻3(Ω), (3) where 𝑢⁰ is the homogenization solution of (1), and dist(𝑥,𝜕Ω) is the distance between the point 𝑥 and the boundary 𝜕Ω.

Unfortunately, the needed CPU time of the method presented in [9] is 𝑂(𝜀^1−𝑛ℎ^−𝑛). In this paper, a high-effective finite element method to compute boundary corrector by virtue of the recovery technique is proposed, and some superconvergence results for the multiscale finite element approximation of (1) are obtained. The rest of this paper is organized as follows. In the next section, we present a multiscale finite element method to compute 𝑢^𝜀(𝑥). Its convergence analysis are shown in Section 3. Finally, some numerical results conforming our analytical estimates are given in Section 4.

Notation. Before closing this section, we would like to fix some notations. First, the Einstein summation is used. Let 𝑄={𝜉∣0<𝜉_𝑖<1,𝑖=1,2}, and the capital letter 𝐶 (with or without subscripts) denotes a positive constant, which is independent of the small parameter 𝜀 and the mesh size ℎ (with or without subscripts).

2. An Improved Multiscale Finite Element Method

Firstly, let us simply recall the homogenization method described in [5].

2.1. Homogenization Method

Let 𝑁_𝑘(𝜉)(𝑘=1,2) be a 1-periodic function, which satisfies 𝜕𝜕𝜉_i𝑎_𝑖𝑗(𝜉)𝜕𝑁_𝑘𝜕𝜉_𝑗𝜕=−𝜕𝜉_𝑖𝑎_𝑖𝑘,inℜ²,_𝑄𝑁_𝑘(𝜉)𝑑𝜉=0. (4) Then, the matrix 𝑎=(𝑎_𝑖𝑗)_2×2 can be obtained by 𝑎_𝑖𝑗=_𝑄𝑎_𝑖𝑗+𝑎_𝑖𝑘𝜕𝑁_𝑗𝜕𝜉_𝑘𝑑𝜉. (5)

The first-order approximation of 𝑢^𝜀(𝑥) can be written as 𝑢(𝑥)=𝑢⁰(𝑥)+𝜀𝑁_𝑘(𝜉)𝜕𝑢⁰(𝑥)𝜕𝑥_𝑘, (6) where 𝑢⁰(𝑥) satisfies the homogenization problem 𝐿₀𝑢⁰(𝑥)𝜕∶=𝜕𝑥_𝑖𝑎_𝑖𝑗𝜕𝑢⁰𝜕𝑥_𝑗(𝑥)𝑢=𝑓,inΩ,₀(𝑥)=𝑔(𝑥),on𝜕Ω. (7)

The boundary corrector term of the homogenization method 𝜃_𝜀 is defined by 𝐿_𝜀𝜃_𝜀𝜃=0,inΩ,_𝜀=−𝜀𝑁_𝑘𝜕𝑢⁰𝜕𝑥_𝑘,on𝜕Ω. (8)

In the next two subsections, we will compute numerically the first-order approximation 𝑢 and the boundary corrector term 𝜃_𝜀, respectively, and furthermore give the multiscale finite element solution of (1).

2.2. Finite Element Approximation of 𝑢

Let 𝒯_ℎ0 be a quasiuniform triangular partition of 𝑄 with the mesh size ℎ₀. 𝑆^ℎ0 denotes the conforming 𝑃₁ finite element spaces with respect to 𝒯_ℎ0, and 𝑆^ℎ0₀=𝑆^ℎ0∩𝐻¹₀(𝑄). The finite element scheme of (4) is to find 𝑁^ℎ0_𝑘∈𝑆^ℎ0 such that _𝑄𝑎_𝑖𝑗(𝜉)𝜕𝑁^ℎ0_𝑘(𝜉)𝜕𝜉_𝑗𝜕𝑣(𝜉)𝜕𝜉_𝑖𝑑𝜉=−_𝑄𝑎_𝑖𝑘(𝜉)𝜕𝑣(𝜉)𝜕𝜉_𝑖𝑑𝜉,∀𝑣∈𝑆^ℎ0₀(𝑄),𝑁^ℎ0_𝑘(𝜉)isa1-periodicfunction. (9) Then, the numerical approximation 𝑎^ℎ0_𝑖𝑗 of 𝑎_𝑖𝑗 can be calculated by 𝑎^ℎ0_𝑖𝑗=_𝑄𝑎_𝑖𝑗(𝜉)+𝑎_𝑖𝑘(𝜉)𝜕𝑁^ℎ0_𝑗(𝜉)𝜕𝜉_𝑘𝑑𝜉. (10)

Let 𝒯_ℎ1 be a quasiuniform triangular partition of Ω with the mesh size ℎ₁ and satisfy min_{𝑒𝑘∈𝒯ℎ1}𝑆_𝑒𝑘≥𝐶ℎ²₁, (11) where 𝑆_𝑒𝑘 is the area of the triangular element 𝑒_𝑘. 𝑆^ℎ1_𝑔 denotes the corresponding conforming 𝑃₁ finite element spaces, and 𝑆^ℎ1₀=𝑆^ℎ1∩𝐻¹₀(𝑄).

The finite element approximation 𝑢^ℎ0,ℎ1₀ of the homogenization problem (7) is to find 𝑢^ℎ0,ℎ1₀∈𝑆^ℎ1_𝑔 such that _Ω𝑎^ℎ0_𝑖𝑗𝜕𝑢^ℎ0,ℎ1₀𝜕𝑥_𝑖𝜕𝑣_ℎ𝜕𝑥_𝑗𝑑𝑥=_Ω𝑓(𝑥)𝑣_ℎ(𝑥)𝑑𝑥,∀𝑣_ℎ∈𝑆^ℎ1₀. (12)

Furthermore, we turn to the computation of 𝜕𝑢^ℎ0,ℎ1₀(𝑥)/𝜕𝑥_𝑘 and 𝑢(𝑥). Let Σ^ℎ1 be the set of all nodal points of the mesh 𝒯_ℎ1. Define 𝑢^ℎ0,ℎ1_𝑘(𝑘=1,2) by the following: (𝐴₁) for all 𝑥∈Σ^ℎ1, the value of 𝑢^ℎ0,ℎ1_𝑘 at the nodal point 𝑥 is the average of𝜕𝑢^ℎ0,ℎ1₀(𝑥)/𝜕𝑥_𝑘 in all elements including 𝑥,(𝐴₂)𝑢^ℎ0,ℎ1_𝑘(𝑥) is a piecewise linear function in every element.

Therefore, we have a numerical approximation 𝑢^ℎ0,ℎ1(𝑥) of 𝑢(𝑥) which is defined by 𝑢^ℎ0,ℎ1(𝑥)=𝑢^ℎ0,ℎ1₀(𝑥)+𝜀𝑁^ℎ0_𝑘(𝜉)𝑢^ℎ0,ℎ1_𝑘(𝑥). (13)

2.3. Finite Element Approximation of 𝜃_𝜀(𝑥)

Let 𝑚 be a positive integer satisfying 2^𝑚−1𝜀<max_𝑥∈Ωdist(𝑥,𝜕Ω)≤2^𝑚𝜀. (14) Then, the domain Ω can be divided by Ω=⋃^𝑚_𝑖=0Ω_𝑖 and Ω_𝑖={()}Ω𝑥∈Ω∣dist𝑥,𝜕Ω≤𝜀,if𝑖=0,_𝑖=𝑥∈Ω∣2^𝑖−1𝜀<dist(𝑥,𝜕Ω)≤2^𝑖𝜀,if1≤𝑖≤𝑚. (15)

Let 𝒯_ℎ be the regular triangular partition of Ω with the mesh size ℎ and satisfy 𝐵₁2^𝑖/2𝜆𝜀ℎ≤𝑆_𝑒𝑖≤2^𝑖/2𝜆⁻¹𝜀ℎ,∀𝑒_𝑖∈𝒯_ℎ,𝑒_𝑖⊂Ω_𝑖∪Ω_𝑖+1,𝐵₂𝑙_𝑖−𝑙_𝑗≤𝐶ℎ𝑙_𝑖,∀𝑒_𝑖∈𝒯_ℎ, (16)(17) where 𝜆 and 𝐶 are independent of 𝑖 and 𝑒_𝑖, 𝑆_𝑒𝑖 denotes the area of 𝑒_𝑖, and 𝑙_𝑖,𝑙_𝑗 are the length of two edges of 𝑒_𝑖. Let Σ^ℎ be the set of all nodal points in 𝒯_ℎ, Σ^ℎ_𝐵=Σ^ℎ∩𝜕Ω, and let 𝑆^ℎ be the conforming 𝑃₁ finite element spaces with respect to 𝒯_ℎ; we define 𝑆^ℎ₁=𝑣∈𝑆^ℎ∣𝑣(𝑥)=−𝜀𝑁^ℎ0_𝑘𝑥𝜀𝑢^ℎ0,ℎ1_𝑘(𝑥),∀𝑥∈Σ^ℎ_𝐵,𝑆^ℎ₀=𝑣∈𝑆^ℎ∣𝑣(𝑥)=0,∀𝑥∈Σ^ℎ_𝐵. (18)

Then, the finite element approximation of 𝜃_𝜀 is to find 𝜃^ℎ0,ℎ1_𝜀^,ℎ∈𝑆^ℎ₁ such that _Ω𝑎_𝑖𝑗(𝜉)𝜕𝜃^ℎ0,ℎ1_𝜀^,ℎ𝜕𝑥_𝑖𝜕𝑣𝜕𝑥_𝑗𝑑𝑥=0,∀𝑣∈𝑆^ℎ₀. (19)

2.4. Multiscale Finite Element Approximation of 𝑢^𝜀

For any 𝑣∈𝑆^ℎ(Ω), we define the linear operator 𝑅_ℎ by 𝜕𝑅_ℎ𝑣(𝑥)𝜕𝑥_𝑖isapiecewisefunctionon𝒯_ℎ,𝜕𝑅_ℎ𝑣(𝑥)𝜕𝑥_𝑖=𝜕𝑣(𝑥)𝜕𝑥_𝑖,if𝑥isthemiddlepointofelementsof𝒯_ℎ. (20) Then, the multiscale finite element approximation 𝑢^{ℎ0,ℎ1,ℎ}(𝑥) of 𝑢^𝜀(𝑥) can be defined by 𝑢^{ℎ0,ℎ1,ℎ}(𝑥)=𝑢^ℎ0,ℎ1(𝑥)+𝑅_ℎ𝜃^ℎ0,ℎ1_𝜀^,ℎ(𝑥). (21)

3. Superconvergence Result of 𝑢^{ℎ0,ℎ1,ℎ}

Firstly, we have the following assumption.

Assumption C1. The functions 𝑎_𝑖𝑗∈𝑊^1,∞(𝑄)∩𝐻²(𝑄), and the homogenization solution 𝑢⁰∈𝑊^2,∞(Ω)∩𝐻⁴(Ω).

Then, we introduce the following lemma.

Lemma 2 (see [1, 5]). Let Ω_𝑟={𝑥∈Ω∣dist(𝑥,𝜕Ω)≥𝑟}. Assuming that (C1) holds, then there exists 𝐶, which is independent of 𝜀 and 𝑟 such that ∇𝜃_{𝜀𝐿2(Ω𝑟)}≤𝐶𝑟^−1/2𝜀𝑢⁰_𝐻2(Ω),𝜃_{𝜀𝐻2(Ω)}≤𝐶𝜀^−1/2𝑢⁰_𝐻3(Ω),𝜃_{𝜀𝐻1(Ω𝑟)}≤𝐶𝑟^−1/2𝜀𝑢⁰_{𝑊1,∞(𝜕Ω)},𝜃_{𝜀𝐻2(Ω𝑟)}≤𝐶𝑟^−1/2𝑢⁰_{𝑊1,∞(𝜕Ω)}. (22)

From Lemma 2, one can easily deduce.

Lemma 3. Assuming that (C1) holds, then there exists 𝐶, which is independent of 𝜀 such that 𝜃_{𝜀𝐻3(Ω)}≤𝐶𝜀^−3/2𝑢⁰_𝐻4(Ω),𝜃_{𝜀𝐻3(Ω𝑟)}≤𝐶𝑟^−1/2𝜀⁻¹𝑢⁰_{𝑊2,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω). (23)(24)

Proof. For 𝑘=1,2, we define 𝑣^𝜀_𝑘(𝑥)=𝜕𝜃_𝜀(𝑥)𝜕𝑥_𝑘. (25) Then, the upper bound of ‖𝑣^𝜀_𝑘‖_𝐻2(Ω) and ‖𝑣^𝜀_𝑘‖_{𝐻1(Ω𝑟)} can be estimated, respectively.
Using the result from (8) and 𝜉=𝑥/𝜀, we have 𝐿_𝜀𝑣^𝜀_𝑘=−𝜀⁻¹𝜕𝜕𝑥_𝑖𝜕𝑎_𝑖𝑗(𝜉)𝜕𝜉_𝑘𝜕𝜃_𝜀𝜕𝑥_𝑗𝑣,𝑥∈Ω,_𝜀𝑘(𝑥)𝜕𝑁=−𝜀^𝑙(𝜉)𝜕𝑢⁰(𝑥)/𝜕𝑥_𝑙𝜕𝑥_𝑘,𝑥∈𝜕Ω. (26) Then, we estimate ‖𝑣^𝜀_𝑘‖_𝐻2(Ω). Obviously, 𝑣^𝜀_𝑘 can be divided into 𝑣^𝜀_𝑘=𝑣^𝜀_𝑘,1+𝑣^𝜀_𝑘,2, (27) where 𝑣^𝜀_𝑘,1 satisfies 𝐿_𝜀𝑣^𝜀_𝑘,1=−𝜀⁻¹𝜕𝜕𝑥_𝑖𝜕𝑎_𝑖𝑗(𝜉)𝜕𝜉_𝑘𝜕𝜃_𝜀𝜕𝑥_𝑗𝑣,𝑥∈Ω,_𝜀𝑘,1=0,𝑥∈𝜕Ω, (28) and 𝑣^𝜀_𝑘,2 satisfies 𝐿_𝜀𝑣^𝜀_𝑘,2𝑣=0,𝑥∈Ω,_𝜀𝑘,2𝜕𝑁=−𝜀^𝑙(𝑥/𝜀)𝜕𝑢⁰/𝜕𝑥_𝑙𝜕𝑥_𝑘,𝑥∈𝜕Ω. (29) Using the result from Lemma 2 and (28), we have 𝑣^𝜀_{𝑘,1𝐻2(Ω)}≤𝐶𝜀⁻¹𝜃_{𝜀𝐻2(Ω)}≤𝐶𝜀^−3/2𝑢⁰_𝐻3(Ω). (30)
Let 𝑣^𝜀_{𝑘,2,𝑘′}=𝜕𝑣^𝜀_𝑘,2/𝜕𝑥_𝑘′, and using the result from (29), we have 𝐿_𝜀𝑣^𝜀_{𝑘,2,𝑘′}𝑣=0,𝑥∈Ω,_{𝜀𝑘,2,𝑘′}𝜕=−𝜀²𝑁_𝑙(𝑥/𝜀)𝜕𝑢⁰/𝜕𝑥_𝑙𝜕𝑥_𝑘𝜕𝑥_𝑘′,𝑥∈𝜕Ω. (31) Following the same line of [5] (1992, Theorem  1.2, pages 124–128), we have 𝑣_{𝑘,2,𝑘′𝐻1(Ω)}≤𝐶𝜀^−3/2𝑢⁰_𝐻4(Ω)≤𝐶𝜀^−3/2𝑢⁰_𝐻4(Ω) (32) which indicates 𝑣^𝜀_{𝑘,2𝐻2(Ω)}≤𝐶𝜀^−3/2𝑢⁰_𝐻4(Ω). (33) Combining (30) with (33), we can derive (24) immediately.
Considering the proof of (23), 𝑣^𝜀_𝑘 can be divided into 𝑣^𝜀_𝑘=𝑣^𝜀_𝑘,1+𝑣^𝜀_𝑘,2, (34) where 𝑣^𝜀_𝑘,1 satisfies 𝐿_𝜀𝑣^𝜀_𝑘,1=−𝜀⁻¹𝜕𝜕𝑥_𝑖𝜕𝑎_𝑖𝑗(𝜉)𝜕𝜉_𝑘𝜕𝜃_𝜀𝜕𝑥_𝑗,𝑥∈Ω_𝑟,𝑣^𝜀_𝑘,1=0,𝑥∈𝜕Ω_𝑟, (35) and 𝑣^𝜀_𝑘,2 satisfies 𝐿_𝜀𝑣^𝜀_𝑘,2=0,𝑥∈Ω_𝑟,𝑣^𝜀_𝑘,2=𝜕𝜃_𝜀𝜕𝑥_𝑘,𝑥∈𝜕Ω_𝑟. (36)
In view of Lemma 2, we have 𝑣^𝜀_{𝑘,1𝐻2(Ω𝑟)}≤𝐶𝜀⁻¹𝜃_{𝜀𝐻2(Ω𝑟)}≤𝐶𝜀⁻¹𝑟^−1/2𝑢⁰_{𝑊2,∞(𝜕Ω)}. (37)
Following the same line of [5] (1992, Theorem  1.2, pages 124–128), we have 𝑣^𝜀_{𝑘,2𝐻2(Ω𝑟)}≤𝐶𝜀⁻¹𝑟^−1/2‖𝜃‖_𝑊^1,∞_(Ω𝑟/2)+𝐶𝑟^−1/2‖𝜃‖_𝑊^2,∞_(Ω𝑟/2)≤𝐶𝜀⁻¹𝑟^−1/2𝑟⁻¹𝜀𝑢⁰_{𝑊2,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω)+𝐶𝑟^−1/2𝑟⁻¹𝑢⁰_{𝑊2,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω)≤𝐶𝜀⁻¹𝑟^−1/2𝑢⁰_{𝑊2,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω). (38) Combining (37) with (38), we can conclude the result of this lemma.

Assuming that 𝒯_ℎ is defined as (16), and let 𝜃^ℎ_𝜀 and 𝜃^𝐼_𝜀 be the linear finite element approximation and the linear interpolation of 𝜃_𝜀 with respect to 𝒯_ℎ, respectively. Then, we have the following.

Lemma 4. Assuming that (C1) holds, then there exists 𝐶 such that ∇𝜃^𝐼_𝜀−𝜃^ℎ_𝜀dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)≤𝐶𝜀^−1/2ℎ²𝑢⁰_{𝑊1,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω). (39)

Proof. Assuming that 𝑚 is defined as (14) and Ω=⋃^𝑚_𝑖=𝑘+1Ω_𝑖 (𝑘<𝑚,𝑘∈ℕ). Considering 𝑎_𝑖𝑗∈𝑊^1,∞(Ω) and using the result from Lemma 2, we have ∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀𝐿2(Ω)}²≤𝐶_Ω𝑎_𝑖𝑗𝑥𝜀𝜕𝜃^ℎ_𝜀,𝑟−𝜃^𝐼_𝜀𝜕𝑥_𝑖𝜕𝜃^ℎ_𝜀−𝜃^𝐼_𝜀𝜕𝑥_𝑗𝑑𝑥≤𝐶_Ω𝑎_𝑖𝑗𝑥𝜀𝜕𝜃_𝜀−𝜃^𝐼_𝜀𝜕𝑥_𝑖𝜕𝜃^ℎ_𝜀,𝑟−𝜃^𝐼_𝜀𝜕𝑥_𝑗≤𝑑𝑥_{𝑘𝑖=1Ω𝑖}𝑎_𝑖𝑗𝑥𝜀𝜕𝜃_𝜀−𝜃^𝐼_𝜀𝜕𝑥_𝑖𝜕𝜃^ℎ_𝜀,𝑟−𝜃^𝐼_𝜀𝜕𝑥_𝑗𝑑𝑥≤𝐶^𝑘_𝑖=12^𝑖/2𝜀ℎ²𝜀⁻¹𝜃_{𝜀𝐻2(Ω𝑖)}∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀,𝑟𝐿2(Ω𝑖)}+𝜃_{𝜀𝐻3(Ω𝑖)}∇𝜃^𝐼_𝜀−𝜃_𝜀^𝑟,ℎ_{𝐿2(Ω𝑖)}≤𝐶^𝑘_𝑖=12^𝑖𝜀²ℎ²𝜀⁻¹(2^𝑖𝜀)^−1/2𝑢⁰_{𝑊1,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω)×∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀,𝑟𝐿2(Ω−Ω𝑟)}≤𝐶2^𝑘/2𝜀^1/2ℎ²𝑢⁰_{𝑊1,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω)×∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀,𝑟𝐿2(Ω)}≤𝐶ℎ²𝑢⁰_{𝑊1,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω)×∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀,𝑟𝐿2(Ω−Ω𝑟)}, (40) which indicates that ∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀𝐿2(Ω)}≤𝐶ℎ²𝑢⁰_{𝑊1,∞(𝜕Ω)}+𝑢⁰_𝐻4(Ω). (41) Then Lemma 4 can be easily derived.

Furthermore, we can obtain the following lemma.

Lemma 5. Assuming that (C1) holds, then there exists 𝐶 such that ∇𝜃^𝐼_𝜀−𝜃^ℎ_𝜀dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)≤𝐶|ln𝜀|ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω). (42)

Proof. Assume that Ω_𝑟𝑗 is defined as Lemma 4. Considering 𝜃^ℎ_𝜀=𝜃^ℎ_𝜀,𝑟0 (∀𝑥∈Ω_𝑟𝑗), and we divide 𝜃^𝐼_𝜀−𝜃^ℎ_𝜀 into 𝜃^𝐼_𝜀−𝜃^ℎ_𝜀=𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀,𝑟𝑗}+^𝑗−1_𝑖=0𝜃^ℎ_{𝜀,𝑟𝑖}−𝜃^ℎ_{𝜀,𝑟𝑖+1}. (43) Then, 𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀𝐻1(Ω−Ω𝑟𝑗)}≤𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀,𝑟𝑗+1𝐻1(Ω−Ω𝑟𝑗)}+^𝑗−1_𝑖=1𝜃^𝐼_{𝜀,𝑟𝑖}−𝜃^ℎ_{𝜀,𝑟𝑖+1𝐻1(Ω−Ω𝑟𝑗)}≤𝑐𝑟_𝑗^1/2ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω)+^𝑗−1_𝑖=1𝜃^𝐼_{𝜀,𝑟𝑖}−𝜃^ℎ_{𝜀,𝑟𝑖+1𝐻1(Ω−Ω𝑟𝑗)}≤𝐶𝑟_𝑗^1/2ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω)+𝐶𝑟_𝑗𝑗ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω)≤𝐶𝑟_𝑗^1/2ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω). (44) Using the result from (44), we have ∇𝜃^𝐼_𝜀−𝜃^ℎ_𝜀dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)≤𝐶^𝑚_𝑗=1∇𝜃^𝐼_𝜀−𝜃^ℎ_𝜀dist(⋅,𝜕Ω)+𝜀_{𝐿2(Ω𝑟𝑗−Ω𝑟𝑗−1)}≤^𝑚_𝑗=1𝐶𝑟_𝑗^−1/2∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀𝐿2(Ω𝑟𝑗−Ω𝑟𝑗−1)}≤^𝑚_𝑗=1𝐶𝑟_𝑗^−1/2∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀𝐿2(Ω−Ω𝑟𝑗−1)}≤^𝑚_𝑗=1𝐶𝑟_𝑗^−1/2𝑟_𝑗^1/2ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω)≤𝐶𝑚ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω)≤𝐶|ln𝜀|ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω). (45)

Based on the previous lemmas, the estimate of ‖|∇(𝜃^𝐼_𝜀−𝜃^ℎ0,ℎ1_𝜀^,ℎ)|/dist(⋅,𝜕Ω)+𝜀‖_𝐿2(Ω) can be given as follows.

Lemma 6. Assuming that (C1) holds, then there exists 𝐶 such that ∇𝜃^𝐼_𝜀−𝜃^ℎ0,ℎ1_𝜀^,ℎdist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)≤𝐶ℎ₀+ℎ₁+ℎ²||×𝑢ln𝜀_{0𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω). (46)

Proof. Assuming that 𝜃^ℎ0,ℎ1_𝜀 satisfies 𝐿_𝜀𝜃^ℎ0,ℎ1_𝜀(𝑥)𝜃=0,𝑥∈Ω,_{ℎ0,ℎ1𝜀}(𝑥)=−𝜀𝑁^ℎ0_𝑘𝑥𝜀𝑢^ℎ0,ℎ1_𝑘(𝑥),𝑥∈𝜕Ω, (47) and 𝜃^{𝐼,ℎ0,ℎ1}_𝜀 is the linear interpolation of 𝜃^ℎ0,ℎ1_𝜀 on 𝒯_ℎ.
Then, we divide 𝜃^𝐼_𝜀−𝜃^ℎ0,ℎ1_𝜀^,ℎ into 𝜃^𝐼_𝜀−𝜃^ℎ0,ℎ1_𝜀^,ℎ=𝜃^𝐼_𝜀−𝜃^{𝐼,ℎ0,ℎ1}_𝜀+𝜃^{𝐼,ℎ0,ℎ1}_𝜀−𝜃^ℎ0,ℎ1_𝜀^,ℎ. (48)
Firstly, considering the first item of the right-hand side of (48) and assuming that 𝜃^ℎ0,ℎ1_𝜀(𝑥) satisfies the problem 𝐿_𝜀𝜃^ℎ0,ℎ1_𝜀(𝑥)𝜃=0,𝑥∈Ω,_{ℎ0,ℎ1𝜀}(𝑥)=−𝜀𝑁^ℎ0_𝑘𝑥𝜀𝑢^ℎ0,ℎ1_𝑘(𝑥),𝑥∈𝜕Ω, (49) we have ∇𝜃^𝐼_𝜀−𝜃^{𝐼,ℎ0,ℎ1}_{𝜀𝐿2(Ω)}≤𝐶∇𝜃_𝜀−𝜃^ℎ0,ℎ1_{𝜀𝐿2(Ω)}≤𝐶∇𝜃_𝜀−𝜃^ℎ0,ℎ1_{𝜀𝐿2(Ω)}+∇𝜃^ℎ0,ℎ1_𝜀−𝜃^ℎ0,ℎ1_{𝜀𝐿2(Ω)}. (50)
Using the same method of Lemma 2, we have 𝜃_𝜀−𝜃^ℎ0,ℎ1_{𝜀𝐿∞(Ω)}≤𝐶𝜀ℎ₀+ℎ₁,∇𝜃_𝜀−𝜃^ℎ0,ℎ1_{𝜀𝐿2(Ω)}≤𝐶𝜀^1/2ℎ₀+ℎ₁. (51)(52)
Then, we have 𝜀𝑁^ℎ0_𝑘𝑥𝜀𝑢^ℎ0,ℎ1_𝑘(𝑥)−𝜀𝑁_𝑘𝑥𝜀𝜕𝑢⁰𝜕𝑥_{𝑘𝐿∞(Ω)}≤𝐶𝜀ℎ₀+ℎ₁𝑢⁰_{𝑊2,∞(Ω)},∇𝜃^ℎ0,ℎ1_𝜀−𝜃^ℎ0,ℎ1_{𝜀𝐿2(Ω)}≤𝐶𝜀^1/2ℎ₀+ℎ₁𝑢⁰_{𝑊2,∞(Ω)}. (53)(54)
Combining (52) with (54), we have ∇𝜃^𝐼_𝜀−𝜃^{𝐼,ℎ0,ℎ1}_{𝜀𝐿2(Ω)}≤𝐶𝜀^1/2ℎ₀+ℎ₁𝑢⁰_{𝑊2,∞(Ω)}. (55)
Next, considering the second item of the right-hand side of (48), 𝜃^{𝐼,ℎ0,ℎ1}_𝜀−𝜃^ℎ0,ℎ1_𝜀^,ℎ can be divided into 𝜃^{𝐼,ℎ0,ℎ1}_𝜀−𝜃^{ℎ,ℎ0,ℎ1}_𝜀=𝜃^{𝐼,ℎ0,ℎ1}_𝜀−𝜃^𝐼_𝜀+𝜃^𝐼_𝜀−𝜃^ℎ_𝜀+𝜃^ℎ_𝜀−𝜃^ℎ0,ℎ1_𝜀^,ℎ. (56)
Similarly, we have ∇𝜃^{𝐼,ℎ0,ℎ1}_𝜀−𝜃^𝐼_{𝜀𝐿2(Ω)}≤𝐶𝜀^1/2ℎ₀+ℎ₁𝑢⁰_{𝑊2,∞(Ω)},∇𝜃^𝐼_𝜀−𝜃^ℎ_{𝜀𝐿2(Ω)}≤𝐶𝜀^1/2ℎ₀+ℎ₁𝑢⁰_{𝑊2,∞(Ω)},∇𝜃^ℎ_𝜀−𝜃^ℎ0,ℎ1_𝜀^,ℎ_𝐿2(Ω)≤𝐶𝜀^1/2ℎ₀+ℎ₁𝑢⁰_{𝑊2,∞(Ω)}. (57)(58) Combining Lemma 6 with (55)–(58), we have (46).

Next, using the extrapolation technique [10], we are in a position to estimate 𝜃_𝜀(𝑥)−𝑅_ℎ4𝜃^ℎ0,ℎ1_𝜀^,ℎ/2(𝑥)−𝜃^ℎ0,ℎ1_𝜀^,ℎ(𝑥)3 (59) instead of 𝜃_𝜀(𝑥)−𝑅_ℎ𝜃^ℎ0,ℎ1_𝜀^,ℎ(𝑥). (60)

Lemma 7. Assuming that (C1) holds, then there exists 𝐶 such that ∇𝜃_𝜀−𝑅_ℎ4𝜃^ℎ0,ℎ1_𝜀^,ℎ/2−𝜃^ℎ0,ℎ1_𝜀^,ℎ/3dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)≤𝐶ℎ₀+ℎ₁+ℎ²|ln𝜀|𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω). (61)

Proof. Let 𝒯_ℎ and 𝜃^𝐼_𝜀 be defined as above. Assuming that 𝜃_𝜀^𝐼,ℎ/2 is the linear interpolation of 𝜃_𝜀 on 𝒯_ℎ/2, we have ∇𝜃_𝜀−𝑅_ℎ4𝜃^ℎ0,ℎ1_𝜀^,ℎ/2−𝜃^ℎ0,ℎ1_𝜀^,ℎ3_𝐿2(Ω)≤∇𝜃_𝜀−𝑅_ℎ4𝜃_𝜀^𝐼,ℎ/2−𝜃^𝐼_𝜀3_𝐿2(Ω)+∇𝑅_ℎ4𝜃_𝜀^𝐼,ℎ/2−𝜃^𝐼_𝜀3−4𝜃^ℎ0,ℎ1_𝜀^,ℎ/2−𝜃^ℎ0,ℎ1_𝜀^,ℎ3_𝐿2(Ω)≤𝐶𝜀^1/2ℎ²|ln𝜀|𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω)+𝐶𝜀^1/2ℎ₀+ℎ₁+ℎ²𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω)≤𝐶𝜀^1/2ℎ₀+ℎ₁+ℎ²|ln𝜀|𝑢⁰_{𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω). (62) Then, (61) can be easily derived.

Next, we turn to estimate ‖∇𝑤^𝜀/dist(⋅,𝜕Ω)+𝜀‖_𝐿²_(Ω).

Lemma 8. Assuming that (C1) holds, then there exists 𝐶 such that ∇𝑤^𝜀dist(⋅,𝜕Ω)+𝜀_𝐿²_(Ω)≤𝐶𝜀|ln𝜀|𝑢⁰_{𝑊2,∞(Ω)}. (63)

Proof. Following the same line of [5] and 𝜔^𝜀=𝑢^𝜀−𝑢−𝜃_𝜀, there exists 𝐶 such that ∇𝑤^𝜀dist(⋅,𝜕Ω)+𝜀_𝐿²_(Ω)≤𝐶∇𝑤^𝜀_{𝐿2|ln𝜀|(Ω)}×1dist(⋅,𝜕Ω)+𝜀_𝐿^{4|ln𝜀|/(2|ln𝜀|−1)}_(Ω)≤𝐶𝜀|ln𝜀|^1/2𝑢⁰_{𝑊2,2|ln𝜀|(Ω)}|ln𝜀|^1/2≤𝐶𝜀|ln𝜀|𝑢⁰_{𝑊2,∞(Ω)}. (64) Finally, noting the definitions of 𝑢, 𝑢^ℎ0,ℎ1, 𝜃_𝜀, 𝑅_ℎ𝜃^ℎ0,ℎ1_𝜀^,ℎ, and 𝜔^𝜀, and combining Lemma 3 with Lemmas 7-8 and Lemma  2.4 in [9], we have ∇𝑢^𝜀−𝑢^{ℎ0,ℎ1,ℎ}dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)≤∇𝑢−𝑢^ℎ0,ℎ1dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)+∇𝜃_𝜀−𝑅_ℎ𝜃^ℎ0,ℎ1_𝜀^,ℎdist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)+∇𝜔^𝜀dist(⋅,𝜕Ω)+𝜀_𝐿²_(Ω)≤𝐶ℎ₁+ℎ₀+𝜀+ℎ²||×𝑢ln𝜀₀‖_𝑊^2,∞_(Ω)+𝑢⁰‖_𝐻⁴_(Ω). (65)

Combining the above lemmas, we can conclude the following result.

Theorem 9. Assuming that (C1) holds, then there exists 𝐶 such that ∇𝑢^𝜀−𝑢^{ℎ0,ℎ1,ℎ}dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)≤𝐶ℎ₁+ℎ₀+𝜀+ℎ²×𝑢|ln𝜀|_{0𝑊2,∞(Ω)}+𝑢⁰_𝐻4(Ω). (66)

4. Numerical Example

In this section, some numerical results will be shown. In order to show the numerical accuracy of the method presented in this paper, the exact solution of problem (1) should firstly be obtained. However, it is very difficult to find them out. Then, the exact solution will be replaced by the finite element solution in a fine mesh with the mesh size 1/256.

It should not be confused that 𝑢^∗ denotes the finite element solution of (1) in a fine mesh, and 𝑢^{ℎ0,ℎ1,ℎ}, obtained by the multiscale finite element scheme presented in the above section, is the multiscale finite element solution of problem (1). Some numerical results will be presented by solving the following model problem: 𝑎₁₁=0.3+2𝑥₁1−𝑥₁,𝑎₁₂=𝑎₂₁=𝑥₁1−𝑥₁𝑥₂1−𝑥₂,𝑎₂₂=0.1+2𝑥₂1−𝑥₂,𝑓(𝑥)=𝑒^𝑥1+𝑥2,𝑔(𝑥)=2sin𝑥₁+4cos𝑥₂,Ω=𝑥∣𝑥₁−0.5²+𝑥₁−0.5²<0.25. (67) Moreover, let 𝑒₀=∇𝑢^∗−𝑢^ℎ0,ℎ1/dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)∇𝑢^∗/dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω),𝑒₁=∇𝑢^∗−𝑢^{ℎ0,ℎ1,ℎ}/dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω)∇𝑢^∗/dist(⋅,𝜕Ω)+𝜀_𝐿2(Ω). (68) In Table 1, the numerical results of the multiscale method for 𝑒₀ and 𝑒₁ are given. It can be seen that the improvement obtained in the final approximation by considering the numerical approximation for the boundary corrector, and the numerical result agree well with the theoretical result from Theorem 9.

Table 1 

Comparison of computational results with 𝜀=1/30.

According to Table 2, it can be seen that ∇𝑢^𝜀(𝑥) can effectively be computed for problem (1) by using the above method, even if dist(𝑥,Ω) is very small. If we only need to get a good numerical solution for problem (1) in Sobolev space 𝐻₁(Ω), the boundary corrector needs not to be computed. However, the boundary corrector is a very important part of error estimate in the real applications. It can be concluded that this method is an exceedingly important and effective finite element algorithm.

Table 2 

Comparison of error order.

Acknowledgments

The authors would like to thank the referees for their valuable suggestions and corrections, which contribute significantly to the improvement of the paper. This research was financially supported by the National Basic Research Program of China (973 Program: 2011CB013800), and the National Natural Science Foundation of China, with Grant nos. 11126132, 50838004, 50908167, and 11101311. The financial support of research is gratefully acknowledged.