Finite Element Approximation of Elliptic Homogenization Problems in Nondivergence-Form

TL;DR

This paper develops and analyzes a finite element numerical scheme for solving elliptic homogenization problems in nondivergence form, providing error estimates and demonstrating effectiveness through numerical experiments.

Contribution

It introduces a new finite element approach for nondivergence-form homogenization problems with rigorous analysis and numerical validation.

Findings

The scheme achieves accurate approximation of solutions.

Numerical experiments confirm the scheme's effectiveness.

Error estimates are established for the proposed method.

Abstract

We use uniform $W^{2,p}$ estimates to obtain corrector results for periodic homogenization problems of the form $A(x/\varepsilon):D^2 u_{\varepsilon} = f$ subject to a homogeneous Dirichlet boundary condition. We propose and rigorously analyze a numerical scheme based on finite element approximations for such nondivergence-form homogenization problems. The second part of the paper focuses on the approximation of the corrector and numerical homogenization for the case of nonuniformly oscillating coefficients. Numerical experiments demonstrate the performance of the scheme.

Equations

$$A\left(\frac{\cdot }{\epsilon }\right):D^2{u}_{\epsilon }:=\sum _{i,j=1}^n{a}_{ij}\left(\frac{\cdot }{\epsilon }\right){\partial }_{ij}^2{u}_{\epsilon }$$

$$u_{\epsilon }=0\text{on }\partial \Omega .$$

$$A\in W^{1,q}(Y)\text{ for some }q>n,$$

$$\{\begin{array}{rl}{A}^0:D^2{u}_0& =f\text{in }\Omega ,\\ \text{hfill}{u}_0& =0\text{on }\partial \Omega .\end{array}$$

$$A^0={\int }_{Y}Am,$$

$$\{\begin{array}{l}{D}^2:(Am)=0\text{in }Y,\\ m\text{is Y-periodic},{\int }_{Y}m=1;\end{array}$$

$$\Vert u_{\epsilon }{\Vert }_{H^{2+s}(\Omega )}=\mathcal{O}\left({\epsilon }^{-s}\right).$$

$$\nabla \cdot \left(A\left(\frac{\cdot }{\epsilon }\right)\nabla u_{\epsilon }\right)+b\left(\frac{\cdot }{\epsilon }\right)\cdot \nabla u_{\epsilon }=f\text{in }\Omega$$

$$\nabla \cdot \left(A\left(\frac{\cdot }{\epsilon }\right)\nabla u_{\epsilon }\right)-\frac{1}{\epsilon }\left(\mathrm{div}A\right)\left(\frac{\cdot }{\epsilon }\right)\cdot \nabla u_{\epsilon }=f\text{in }\Omega ,$$

$$\Vert m-m_h{\Vert }_{L^2(Y)}+h\Vert m-m_h{\Vert }_{H^1(Y)}\lesssim h\underset{{\tilde{v}}_h\in {\tilde{M}}_h}{\inf }\Vert m-({\tilde{v}}_h+1){\Vert }_{H^1(Y)},$$

$$A_h^0:={\int }_Y{\mathcal{I}}_h(A{m}_h),$$

$$\{\begin{array}{rl}{A}_h^0:D^2{u}_0^h& =f\text{in }\Omega ,\\ \text{hfill}{u}_0^h& =0\text{on }\partial \Omega ,\end{array}$$

$${\Vert u_0^h-u_0^{h,k}\Vert }_{H^s(\Omega )}\lesssim k\Vert f{\Vert }_{H^{s-1}(\Omega )},$$

$${\Vert u_0-u_0^{h,k}\Vert }_{H^s(\Omega )}\lesssim (h+k)\Vert f{\Vert }_{H^{s-1}(\Omega )},$$

$${\Vert u_{\epsilon }-u_0-{\epsilon }^2\sum _{i,j=1}^n{\chi }_{ij}\left(\frac{\cdot }{\epsilon }\right){\partial }_{ij}^2{u}_0\Vert }_{H^2(\Omega )}\lesssim \sqrt{\epsilon }\Vert u_0{\Vert }_{W^{2,\infty }(\Omega )}+\epsilon \Vert u_0{\Vert }_{H^4(\Omega )},$$

$$\{\begin{array}{l}A:D^2{\chi }_{ij}=a_{ij}^0-a_{ij}\text{in }Y,\\ {\chi }_{ij}\text{is Y-periodic},{\int }_Y{\chi }_{ij}=0.\end{array}$$

$$\Vert u_{\epsilon }-u_0{\Vert }_{H^1(\Omega )}+\sum _{k,l=1}^n{\Vert {\partial }_{kl}^2{u}_{\epsilon }-\left({\partial }_{kl}^2{u}_0+\sum _{i,j=1}^n\left({\partial }_{kl}^2{\chi }_{ij}\right)\left(\frac{\cdot }{\epsilon }\right){\partial }_{ij}^2{u}_0\right)\Vert }_{L^2(\Omega )}=\mathcal{O}(\sqrt{\epsilon }),$$

$$\{\begin{array}{rl}A\left(\cdot ,\frac{\cdot }{\epsilon }\right):D^2{u}_{\epsilon }& =f\text{in }\Omega ,\\ \text{hfill}{u}_{\epsilon }& =0\text{on }\partial \Omega ,\end{array}$$

$$A\in W^{2,\infty }(\Omega ;W^{1,q}(Y))\text{ for some }q>n.$$

$${\Vert u_{\epsilon }-u_0-{\epsilon }^2\sum _{i,j=1}^n{\chi }_{ij}\left(\cdot ,\frac{\cdot }{\epsilon }\right){\partial }_{ij}^2{u}_0\Vert }_{H^2(\Omega )}\lesssim \sqrt{\epsilon }\Vert u_0{\Vert }_{W^{2,\infty }(\Omega )}+\epsilon \Vert u_0{\Vert }_{H^4(\Omega )},$$

$$Y:=(0,1{)}^n,$$

$$A=A^{\mathrm{T}}:{\mathbb{R}}^n\to {\mathbb{R}}^{n\times n}$$

$$\{\begin{array}{l}A\in W^{1,q}(Y)\text{ for some }q\in (n,\infty ],\\ A\text{ is }\text{Y-periodic},\\ \exists \lambda ,\Lambda >0:\lambda |\xi {|}^2\le A(y)\xi \cdot \xi \le \Lambda |\xi {|}^2\forall \xi ,y\in {\mathbb{R}}^n.\end{array}$$

$$A\in C^{0,\alpha }({\mathbb{R}}^n)\text{ for some }0<\alpha \le 1.$$

$$\{\begin{array}{rl}A\left(\frac{\cdot }{\epsilon }\right):D^2{u}_{\epsilon }& =f\text{in }\Omega ,\\ \text{hfill}{u}_{\epsilon }& =0\text{on }\partial \Omega ,\end{array}$$

$$(\Omega ,A,f)\in {\mathcal{G}}^{m,p}⟺\{\begin{array}{r}\Omega \subset {\mathbb{R}}^n\text{ is a bounded C2,γ domain},\gamma \in (0,1),\\ A=A^{\mathrm{T}}:{\mathbb{R}}^n\to {\mathbb{R}}^{n\times n}\text{ satisfies eqrefAass,}\\ f\in W^{m,p}(\Omega ),\end{array}$$

$$(\Omega ,A,f)\in {\mathcal{H}}^m⟺\{\begin{array}{r}\Omega \subset {\mathbb{R}}^n\text{ is a bounded convex domain,}\\ A=A^{\mathrm{T}}:{\mathbb{R}}^n\to {\mathbb{R}}^{n\times n}\text{ satisfies eqrefAass,}\\ \exists \delta \in (0,1]:\frac{|A{|}^2}{(\mathrm{tr}A{)}^2}\le \frac{1}{n-1+\delta }\text{ in }{\mathbb{R}}^n,\\ f\in H^m(\Omega ).\end{array}$$

$$\frac{|A(y){|}^2}{(\mathrm{tr}A(y){)}^2}\le \frac{1}{n-1+\delta }\forall y\in {\mathbb{R}}^n,$$

$$\frac{|A(y){|}^2}{(\mathrm{tr}A(y){)}^2}=1-\frac{2\det A(y)}{(\mathrm{tr}A(y){)}^2}\le 1-\frac{2{\lambda }^2}{4{\Lambda }^2}=\frac{1}{1+\delta }\forall y\in {\mathbb{R}}^n$$

$$(\Omega ,A,f)\in {\mathcal{H}}^m⟺\{\begin{array}{r}\Omega \subset {\mathbb{R}}^n\text{ is a bounded convex domain,}\\ A=A^{\mathrm{T}}:{\mathbb{R}}^n\to {\mathbb{R}}^{n\times n}\text{ satisfies eqrefAass,}\\ f\in H^m(\Omega ).\end{array}$$

Full text

1. 1 Introduction

2. 2 Homogenization of Elliptic Problems in Nondivergence-Form

3. 2.1 Framework

4. 2.2 Transformation into Divergence-Form

5. 2.3 Uniform $W^{2,p}$ Estimates and Homogenization Theorem

6. 2.4 Correctors

7. 3 The Numerical Scheme

8. 3.1 Numerical Homogenization Scheme

9. 3.1.1 Approximation of $m$

10. 3.1.2 Approximation of $A^{0}$

11. 3.1.3 Approximation of $u_{0}$

12. 3.2 Approximation of the Corrector

13. 3.3 Approximation of $u_{\varepsilon}$

14. 3.4 Nonuniformly Oscillating Coefficients

15. 4 Numerical Experiments

16. 4.1 Problem with a Known $u_{0}$

17. 4.2 Problem with an Unknown $u_{0}$

18. 4.3 Nonuniformly Oscillating Coefficients

19. 5 Conclusion

1. Introduction

In this work we consider second-order elliptic equations of nondivergence structure, involving rapidly oscillating coefficients, of the form

[TABLE]

subject to the homogeneous Dirichlet boundary condition

[TABLE]

Here we assume that $\Omega\subset\mathbb{R}^{n}$ is a sufficiently regular bounded domain, $\varepsilon>0$ is small, and that $A=(a_{ij}):\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ is a symmetric, uniformly elliptic and $(0,1)^{n}$-periodic matrix-valued function such that

[TABLE]

where $Y:=(0,1)^{n}$ denotes the unit cell, see (2.1). The main goal of this paper is to propose and analyze a numerical homogenization scheme for (1.1), (1.2) that is based on finite element approximations.

The theory of periodic homogenization is concerned with the limiting behavior of the solutions as the oscillation parameter $\varepsilon$ tends to zero. For the problem (1.1), (1.2) under consideration a classical homogenization theorem (see [6, Sec. 3, Theorem 5.2]) states that the solution sequence $(u_{\varepsilon})_{\varepsilon>0}$ converges in an appropriate Sobolev space to the solution $u_{0}$ to the problem

[TABLE]

Here $A^{0}\in\mathbb{R}^{n\times n}$ is the constant matrix given by

[TABLE]

and $m:\mathbb{R}^{n}\rightarrow\mathbb{R}$ is the invariant measure, i.e. the solution to the problem

[TABLE]

see Section 2 for further details. The task of numerical homogenization is the numerical approximation of the matrix $A^{0}$ and the solution $u_{0}$ to the homogenized problem (1.3). As it turns out, $u_{0}$ provides a good approximation to $u_{\varepsilon}$ in $H^{1}(\Omega)$, and by adding corrector terms it is possible to obtain an $H^{2}(\Omega)$-norm approximation. Note that the approximation of (1.1), (1.2) by a standard $H^{2}(\Omega)$-conforming finite element method does not yield error bounds independent of $\varepsilon$, since for $s>0$ one has that

[TABLE]

The motivation for investigating second-order elliptic problems in nondivergence-form comes from physics, engineering, as well as mathematical areas such as stochastic analysis. A notable example of a nonlinear PDE of nondivergence structure is the Hamilton–Jacobi–Bellman equation, which arises in stochastic control theory. The asymptotic behavior of PDEs with rapidly oscillating coefficients is also of importance when micro-inhomogeneous media are investigated.

Over the past decades significant work has been done on periodic homogenization of elliptic problems in divergence-form; numerical homogenization for nondivergence-form problems is however less developed.

The theory of homogenization of divergence-form problems such as

[TABLE]

with periodic and sufficiently regular $A:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ and $b:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n}$ is extensively covered in the books [1, 6, 9, 20]. For divergence-form problems, various multiscale finite element methods (MsFEM) have been developed, which have the advantage over classical finite element methods of providing accurate approximations for very small values of $\varepsilon$ even for moderate values of the grid size. The book [10] by Efendiev and Hou contains a detailed overview of these methods.

It is important to note that although, if $A$ is sufficiently smooth, equation (1.1) can be rewritten in divergence-form,

[TABLE]

this equation does not fit into the framework of divergence-form homogenization problems such as (1.5), because of the $\varepsilon^{-1}$ term in front of the first-order term in (1.6).

For the theory of homogenization of nondivergence-form problems such as (1.1) we refer to the monograph [6] by Bensoussan, Lions and Papanicolaou, to the paper [2] by Avellaneda and Lin, and the references therein. In [5], Bensoussan, Boccardo and Murat study the more general problem involving a Hamiltonian with quadratic growth. Numerical homogenization for nondivergence-form problems using finite difference schemes has been considered in [11] by Froese and Oberman.

The first step in the development of the proposed numerical homogenization scheme is the construction of a finite element method to obtain approximations $(m_{h})_{h>0}\subset H^{1}_{\text{per}}(Y)$ to the invariant measure with optimal order convergence rate

[TABLE]

where $\tilde{M}_{h}$ denotes the finite-dimensional subspace of $H^{1}_{\text{per}}(Y)$ consisting of continuous $Y$-periodic piecewise linear functions on the triangulation with zero mean over $Y$; see Theorem 3.1.

Throughout this work, we use the notation $a\lesssim b$ for $a,b\in\mathbb{R}$ to denote that $a\leq Cb$ for some constant $C>0$ that does not depend on $\varepsilon$ and the discretization parameters.

The second step is to obtain approximations $(A^{0}_{h})_{h>0}\subset\mathbb{R}^{n\times n}$ to the constant matrix $A^{0}$; see Lemma 3.1. To this end, the integrand in (1.4) is replaced by its continuous piecewise linear interpolant and the invariant measure $m$ is replaced by the approximation $m_{h}$, i.e.,

[TABLE]

which can be computed exactly using an appropriate quadrature rule.

The third step is to perform an $H^{s}(\Omega)$-conforming ($s\in\{1,2\}$) finite element approximation for the problem

[TABLE]

on a family of triangulations of the computational domain $\Omega$, parametrized by a discretization parameter $k>0$, measuring the granularity of the triangulation, to obtain $(u_{0}^{h,k})_{h,k>0}\subset H^{s}(\Omega)\cap H^{1}_{0}(\Omega)$ with

[TABLE]

where the constant is independent of $h$; see Lemma 3.3. Note that for the sake of approximating $u_{0}$, an $H^{1}(\Omega)$-conforming finite element method is sufficient.

The approximation $(u_{0}^{h,k})_{h,k>0}\subset H^{s}(\Omega)\cap H^{1}_{0}(\Omega)$ obtained by this procedure approximates $u_{0}$, i.e., the solution to (1.3), with convergence rate

[TABLE]

which can be improved to $\mathcal{O}(h^{2}+k)$ for more regular $A$; see Theorem 3.2, Theorem 3.3 and Remark 3.3.

Concerning the approximation of $u_{\varepsilon}$, i.e., the solution to (1.1), (1.2), we show in Section 2 that under certain assumptions on the domain and the right-hand side, one has that

[TABLE]

where the corrector functions $\chi_{ij}:\mathbb{R}^{n}\rightarrow\mathbb{R}$, $i,j=1,\dots,n$, are defined as the solutions to

[TABLE]

This provides us with the estimate

[TABLE]

which shows that $u_{0}$ is a good $H^{1}(\Omega)$ approximation to $u_{\varepsilon}$ for small $\varepsilon$, and we show in Sections 3.2 and 3.3 how the above estimate can be used to obtain approximations to $D^{2}u_{\varepsilon}$. Note that in order to approximate $u_{\varepsilon}$ in the $H^{1}(\Omega)$-norm, it is sufficient to approximate $u_{0}$ in the $H^{1}(\Omega)$-norm. However, for an approximation of $D^{2}u_{\varepsilon}$ based on the above corrector estimate, we need to approximate $u_{0}$ in the $H^{2}(\Omega)$-norm.

In Section 3.4, we extend our results to the case of nonuniformly oscillating coefficients, i.e., to problems of the form

[TABLE]

where $A=A(x,y):\Omega\times\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ is a symmetric, uniformly elliptic matrix-valued function that is $Y$-periodic in $y$ for fixed $x\in\Omega$, and such that

[TABLE]

We prove the corrector estimate

[TABLE]

where $u_{0}$ is the solution to the homogenized problem corresponding to (1.7) and $\chi_{ij}$ are certain corrector functions. We then discuss the numerical approximation of $u_{\varepsilon}$ based on this corrector estimate, see Section 3.4.

2. Homogenization of Elliptic Problems in Nondivergence-Form

2.1. Framework

We denote the unit cell in $\mathbb{R}^{n}$ by

[TABLE]

and consider a symmetric matrix-valued function

[TABLE]

with the properties

[TABLE]

By Sobolev embedding, we then have that

[TABLE]

For $\varepsilon>0$, we are concerned with the problem

[TABLE]

where the triple $(\Omega,A,f)$ satisfies one of the following sets of assumptions.

Definition 2.1** (Sets of assumptions $\mathcal{G}^{m,p}$, $\mathcal{H}^{m}$).**

For $m\in\mathbb{N}_{0}$ and $p\in(1,\infty)$, we define the set of assumptions $\mathcal{G}^{m,p}$ as

[TABLE]

and the set of assumptions $\mathcal{H}^{m}$ as

[TABLE]

Remark 2.1**.**

For $n=2$, the Cordes condition, i.e., that there exists a $\delta\in(0,1]$ such that

[TABLE]

is a consequence of the uniform ellipticity condition. Indeed, for $A=A^{\mathrm{T}}:\mathbb{R}^{2}\rightarrow\mathbb{R}^{2\times 2}$ satisfying (2.1), we have that

[TABLE]

with $\delta=\frac{\lambda^{2}}{2\Lambda^{2}-\lambda^{2}}\in(0,1]$. Therefore, when $n=2$, the set $\mathcal{H}^{m}$ can be simplified to

[TABLE]

The following theorem asserts well-posedness of the problem (2.2); see [12, Theorem 9.15] and [19, Theorem 3].

Theorem 2.1** (Existence and uniqueness of strong solutions).**

Assume either that $(\Omega,A,f)\in\mathcal{G}^{0,p}$ for some $p\in(1,\infty)$, or that $(\Omega,A,f)\in\mathcal{H}^{0}$ and $p=2$. Then, for any $\varepsilon>0$, the problem (2.2) admits a unique solution $u_{\varepsilon}\in W^{2,p}(\Omega)\cap W^{1,p}_{0}(\Omega)$.

2.2. Transformation into Divergence-Form

We recall a well-known procedure to transform the problem (2.2) into divergence-form; see [2, 6]. We use the notation

[TABLE]

Let us start by introducing the notion of invariant measure; see [6].

Lemma 2.1** (Invariant measure and solvability condition).**

Let $A=A^{\mathrm{T}}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ satisfy (2.1). Then, there exists a unique solution $m:\mathbb{R}^{n}\rightarrow\mathbb{R}$ to the problem

[TABLE]

The function $m$ is called the invariant measure. There holds $m\in W^{1,q}(Y)$, see [7, 8], and there exist constants $\bar{m},M>0$ such that

[TABLE]

Moreover, for a $Y$-periodic function $g\in L^{2}(\mathbb{R}^{n})$, the (adjoint) problem

[TABLE]

admits a solution $u\in W_{\text{per}}(Y)$ if and only if

[TABLE]

With the invariant measure at hand, we can easily convert the problem into divergence-form as follows. We define a matrix-valued function $B=(b_{ij})_{1\leq i,j\leq n}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ by

[TABLE]

with $v_{l}\in W_{\text{per}}(Y)$ denoting the solution to

[TABLE]

for $1\leq l\leq n$. Since $A\in W^{1,q}(Y)$ and $m\in W^{1,q}(Y)$, by elliptic regularity one has that $v_{l}\in W^{2,q}(Y)$ for any $1\leq l\leq n$. Hence, we have

[TABLE]

Further, we observe that $B$ is skew-symmetric, $Y$-periodic with zero mean over $Y$, and that

[TABLE]

Now we let

[TABLE]

Then, since

[TABLE]

and using the fact that $B$ is skew-symmetric, we obtain

[TABLE]

i.e., we have converted (2.2) into divergence-form

[TABLE]

and it is straightforward to check that $A^{\mathrm{div}}$ is $Y$-periodic, Hölder continuous on $\mathbb{R}^{n}$ and uniformly elliptic.

2.3. Uniform $W^{2,p}$ Estimates and Homogenization Theorem

The transformation described in the previous section can be used to obtain uniform $W^{2,p}(\Omega)$ a priori estimates for the solution of (2.2), which are crucial in deriving homogenization results.

Theorem 2.2** (Uniform $W^{2,p}$ a priori estimates).**

Assume either that $(\Omega,A,f)\in\mathcal{G}^{0,p}$ for some $p\in(1,\infty)$, or that $(\Omega,A,f)\in\mathcal{H}^{0}$ and $p=2$. Then, for $\varepsilon\in(0,1]$, the solution $u_{\varepsilon}\in W^{2,p}(\Omega)\cap W^{1,p}_{0}(\Omega)$ to (2.2), whose existence and uniqueness are guaranteed by Theorem 2.1, satisfies

[TABLE]

with the constant absorbed into the notation $\lesssim$ being independent of $\varepsilon$.

Proof.

Let us first assume that $(\Omega,A,f)\in\mathcal{G}^{0,p}$ for some $p\in(1,\infty)$. We showed in the previous section that we can transform problem (2.2) into the divergence-form problem (2.6), where $A^{\mathrm{div}}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ is a $Y$-periodic, Hölder continuous, and uniformly elliptic matrix-valued function satisfying

[TABLE]

Therefore, we can apply [3, Theorem D] to problem (2.6) to obtain

[TABLE]

with constants independent of $\varepsilon$, where we have used the property (2.4) of the invariant measure in the second inequality.

Let us now assume that $(\Omega,A,f)\in\mathcal{H}^{0}$. Noting that (2.3) implies the Cordes condition for $A\left(\frac{\cdot}{\varepsilon}\right)$ with the same constant $\delta\in(0,1]$ for any $\varepsilon>0$, the proof of [19, Theorem 3] yields the estimate

[TABLE]

where $\gamma$ is the function given by

[TABLE]

We observe that by (2.1), there exist constants $\bar{\gamma},\Gamma>0$ such that

[TABLE]

Therefore, we obtain from (2.7) the bound

[TABLE]

with a constant that is independent of $\varepsilon$. ∎

This leads to a simple proof of the homogenization theorem for problem (2.2), using the compactness of the embedding $W^{2,p}(\Omega)\hookrightarrow W^{1,p}(\Omega)$ and the fact that we can rewrite the problem as (2.6).

Theorem 2.3** (Homogenization theorem for nondivergence-form problems).**

Assume either that $(\Omega,A,f)\in\mathcal{G}^{0,p}$ for some $p\in(1,\infty)$, or that $(\Omega,A,f)\in\mathcal{H}^{0}$ and $p=2$. Then the solution $u_{\varepsilon}\in W^{2,p}(\Omega)\cap W^{1,p}_{0}(\Omega)$ to (2.2) converges weakly in $W^{2,p}(\Omega)$ to the solution $u_{0}\in W^{2,p}(\Omega)\cap W^{1,p}_{0}(\Omega)$ of the homogenized problem

[TABLE]

with $A^{0}=(a^{0}_{ij})_{1\leq i,j\leq n}\in\mathbb{R}^{n\times n}$ being the constant matrix whose entries are given by

[TABLE]

where $m$ is the invariant measure, as defined in Lemma 2.1.

Proof.

By Theorem 2.2, the reflexivity of $W^{2,p}(\Omega)$, the compactness of the embedding $W^{2,p}(\Omega)\hookrightarrow W^{1,p}(\Omega)$, and the properties of the trace operator, there exists a $u_{0}\in W^{2,p}(\Omega)\cap W^{1,p}_{0}(\Omega)$ such that (for a subsequence, not indicated,)

[TABLE]

We can transform (2.2) as in Section 2.2 into the divergence-form problem (2.6) with

[TABLE]

being $Y$-periodic, Hölder continuous and uniformly elliptic on $\mathbb{R}^{n}$. Recalling that $B$ is of mean zero over $Y$, we have

[TABLE]

Since we have that

[TABLE]

we can pass to the limit in the weak formulation of (2.6) to obtain that $u_{0}\in W^{2,p}(\Omega)\cap W^{1,p}_{0}(\Omega)$ solves (2.8). We conclude the proof by noting that (2.8) admits a unique strong solution in $W^{2,p}(\Omega)\cap W^{1,p}_{0}(\Omega)$. ∎

2.4. Correctors

We show that by adding corrector terms to the solution $u_{0}$ of the homogenized problem, we obtain a $W^{2,p}$ convergence result.

Theorem 2.4** (Corrector estimate I).**

Assume either that $(\Omega,A,f)\in\mathcal{G}^{2,p}$ for some $p\in(1,\infty)$, or that $(\Omega,A,f)\in\mathcal{H}^{2}$ and $p=2$. Let $\varepsilon\in(0,1]$ and assume that

[TABLE]

Introducing the corrector function $\chi_{ij}$, $1\leq i,j\leq n$, as the solution to

[TABLE]

and a boundary corrector $\theta_{\varepsilon}$, as the solution to

[TABLE]

the following bound holds:

[TABLE]

Proof.

First, we note that since $A\in C^{0,\alpha}(\mathbb{R}^{n})$, we have $\chi_{ij}\in C^{2,\alpha}(\mathbb{R}^{n})$ for any $1\leq i,j\leq n$ by elliptic regularity theory. A direct computation shows that the function

[TABLE]

solves the problem

[TABLE]

where

[TABLE]

Note that since $u_{0}\in W^{4,p}(\Omega)$, one has that

[TABLE]

with the constant being independent of $\varepsilon$. We then have that $d_{\varepsilon}:=\tilde{u}_{\varepsilon}-u_{\varepsilon}$ satisfies

[TABLE]

Therefore, by the definition of the boundary corrector,

[TABLE]

We conclude using the estimate from Theorem 2.2 that

[TABLE]

and (2.10) holds. ∎

The following theorem shows that if $u_{0}\in W^{4,p}(\Omega)\cap W^{2,\infty}(\Omega)$, then we can absorb the term involving the boundary corrector into the right-hand side at the cost of powers of $\varepsilon$.

Theorem 2.5** (Corrector estimate II).**

Assume either that $(\Omega,A,f)\in\mathcal{G}^{2,p}$ for some $p\in(1,\infty)$, or that $(\Omega,A,f)\in\mathcal{H}^{2}$ and $p=2$. Let $\varepsilon\in(0,1]$ and assume that

[TABLE]

Then,

[TABLE]

Proof.

Let $\eta\in C_{c}^{\infty}(\mathbb{R}^{n})$ be a cut-off function with $0\leq\eta\leq 1$,

[TABLE]

and let $\eta$ satisfy

[TABLE]

We introduce the function

[TABLE]

and verify that

[TABLE]

where $S_{1},S_{2}$ and $S_{3}$ are given by

[TABLE]

Therefore, $\tilde{\theta}_{\varepsilon}$ satisfies

[TABLE]

Since $u_{0}\in W^{4,p}(\Omega)\cap W^{2,\infty}(\Omega)$ by assumption, the right-hand side belongs to $L^{p}(\Omega)$, and we have by Theorem 2.2 that

[TABLE]

We look at the terms on the right-hand side separately and start with $S_{1}$. Using the boundedness of $A$ and the fact that $\chi_{ij}\in W^{2,\infty}(\mathbb{R}^{n})$, we have

[TABLE]

For $S_{2}$, we obtain similarly that

[TABLE]

Finally, for $S_{3}$, we have that

[TABLE]

Altogether, we have shown that

[TABLE]

By direct computation, using the bounds

[TABLE]

we can show that

[TABLE]

Therefore, using the triangle inequality, we obtain that

[TABLE]

We conclude that

[TABLE]

The claim now follows from (2.10). ∎

Let us remark that $W^{4,p}(\Omega)\hookrightarrow W^{2,\infty}(\Omega)$ for $p>\frac{n}{2}$, i.e., assumption (2.11) is a consequence of $u_{0}\in W^{4,p}(\Omega)$; in particular, for dimensions $n\in\{2,3\}$ and $p=2$, one can replace condition (2.11) by $u_{0}\in H^{4}(\Omega)$.

Let us recall that $u_{0}$ is the solution to the elliptic constant-coefficient problem (2.8). For bounded convex polygonal domains ($n=2$), $u_{0}\in H^{4}(\Omega)$ can be ensured by assuming that $f\in H^{2}(\Omega)$ satisfies certain compatibility conditions at the corners of the domain. In the case of Poisson’s equation on $\Omega=(0,1)^{2}$, a necessary and sufficient condition for $u_{0}\in H^{4}(\Omega)\cap H^{1}_{0}(\Omega)$ is that $f\in H^{2}(\Omega)$ and $f=0$ at the corners of $\Omega$, see [15]. We note that these conditions are satisfied for functions $f\in H^{2}(\Omega)$ such that $\mathrm{supp}(f)\Subset\Omega$, see [13].

3. The Numerical Scheme

3.1. Numerical Homogenization Scheme

The first step is to approximate the invariant measure.

3.1.1. Approximation of $m$

For the approximation of the invariant measure $m$, we consider a shape-regular triangulation of $Y$ into triangles with longest edge $h>0$ and let

[TABLE]

be the finite-dimensional subspace of $W_{\text{per}}(Y)$ consisting of continuous $Y$-periodic piecewise linear functions on the triangulation with zero mean over $Y$. We assume that

[TABLE]

Then we have the following approximation result for $m$.

Theorem 3.1** (Approximation of the invariant measure).**

Let $A=A^{\mathrm{T}}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ satisfy (2.1). Then, for $h>0$ sufficiently small, there exists a unique $\tilde{m}_{h}\in\tilde{M}_{h}$ such that

[TABLE]

and writing

[TABLE]

we have that

[TABLE]

where $m$ is the invariant measure, as defined in Lemma 2.1.

Remark 3.1**.**

In particular, since

[TABLE]

we have that

[TABLE]

as $h$ tends to zero.

Proof of Theorem 3.1.

We observe that $m=\tilde{m}+1$, where $\tilde{m}$ is the unique solution to the problem

[TABLE]

i.e.,

[TABLE]

where

[TABLE]

We further observe that (3.1) is equivalent to

[TABLE]

We start by showing boundedness of $a$ and a Gårding-type inequality. We claim that there exist constants $C_{b},C_{g}>0$ such that

[TABLE]

and

[TABLE]

Let us first show (3.3). For $u,v\in W_{\text{per}}(Y)$, by Hölder’s inequality and Sobolev embeddings (note that, according to (2.1), $q>n$), we have that

[TABLE]

Using the fact that $A\in W^{1,q}(Y)\hookrightarrow L^{\infty}(Y)$ since $q>n$, we obtain the bound

[TABLE]

for any $u,v\in W_{\text{per}}(Y)$, i.e. (3.3) holds.

Let us now show the estimate (3.4). For $u\in W_{\text{per}}(Y)$, by ellipticity and Hölder’s inequality, we have

[TABLE]

For the second term we use the Gagliardo–Nirenberg inequality and Young’s inequality to obtain

[TABLE]

Therefore, we have

[TABLE]

for any $u\in W_{\text{per}}(Y)$, i.e., (3.4) holds with

[TABLE]

We use Schatz’s method to derive an a priori estimate; see [18].

From our Gårding-type inequality (3.4) we see that (note $\tilde{m}-\tilde{m}_{h}\in W_{\text{per}}(Y)$)

[TABLE]

By Galerkin-orthogonality and boundedness, we have for any $\tilde{v}_{h}\in\tilde{M}_{h}$ that

[TABLE]

and taking the infimum over all $\tilde{v}_{h}\in\tilde{M}_{h}$, we find

[TABLE]

Combining this estimate with (3.5) yields

[TABLE]

Next, we use an Aubin–Nitsche-type duality argument.

Let $\phi\in W_{\text{per}}(Y)$ be the unique solution to

[TABLE]

We note that the solvability condition (2.5) is satisfied:

[TABLE]

We have, using the bounds on the invariant measure (2.4), the weak formulation of (3.7) and the symmetry of $A$, that

[TABLE]

Next, we use Galerkin orthogonality, the boundedness (3.3) and an interpolation inequality to obtain

[TABLE]

where $\mathcal{I}_{h}\phi$ denotes the continuous piecewise linear interpolant of $\phi$ on the triangulation. Finally, by a regularity estimate for $\phi$ and the bounds on the invariant measure (2.4), we arrive at the bound

[TABLE]

which provides us with the estimate

[TABLE]

for some $C_{0}>0$. Combining this with (3.6) we have

[TABLE]

Therefore, for $h$ sufficiently small, we arrive at the bounds

[TABLE]

and

[TABLE]

We have thus established the a priori estimate

[TABLE]

which immediately implies existence and uniqueness of solutions to (3.2).

Finally, using that $m=\tilde{m}+1$ and $m_{h}=\tilde{m}_{h}+1$, we conclude that

[TABLE]

∎

3.1.2. Approximation of $A^{0}$

We use this finite element approximation of the invariant measure to obtain an approximation to the constant matrix

[TABLE]

To this end, we first replace the invariant measure $m$ by the approximation $m_{h}$ from Theorem 3.1, and then replace the integrand by its piecewise linear interpolant,

[TABLE]

This integral can be computed exactly using an appropriate quadrature rule. The following lemma gives an error estimate for this approximation.

Lemma 3.1** (Approximation of $A^{0}$).**

Let $A=A^{\mathrm{T}}:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ satisfy (2.1). Further, let $A^{0}=(a_{ij}^{0})\in\mathbb{R}^{n\times n}$ be the constant matrix given by Theorem 2.3, let $m_{h}$ be the approximation to the invariant measure given by Theorem 3.1, and let $A_{h}^{0}=(a_{ij,h}^{0})\in\mathbb{R}^{n\times n}$ be the matrix given by

[TABLE]

Then, for $h>0$ sufficiently small, $A_{h}^{0}$ is elliptic and

[TABLE]

Proof.

Fix $1\leq i,j\leq n$. Using the definition of $A^{0}=(a_{ij}^{0})$, i.e.,

[TABLE]

we obtain the estimate

[TABLE]

For the first term, we have

[TABLE]

For the second term, let us first note that using $a_{ij}\in W^{1,q}(Y)$ with $q>n$ and Sobolev embeddings, we have

[TABLE]

Therefore, using a standard interpolation error bound, we obtain

[TABLE]

By Theorem 3.1, for $h>0$ sufficiently small, we have that

[TABLE]

Finally, we note that this implies that for $h>0$ sufficiently small, $A_{h}^{0}$ is elliptic. ∎

3.1.3. Approximation of $u_{0}$

For the approximation of the solution $u_{0}$ to the homogenized problem, we use the following comparison result for the error committed when replacing $A^{0}$ by $A^{0}_{h}$.

Lemma 3.2** (Comparison result).**

Assume either that $(\Omega,A,f)\in\mathcal{G}^{0,2}$ or that $(\Omega,A,f)\in\mathcal{H}^{0}$. Let $A_{h}^{0}\in\mathbb{R}^{n\times n}$ be the approximation to $A^{0}$ as in Lemma 3.1. Then, for $h>0$ sufficiently small, we have that

[TABLE]

where $u_{0}^{h}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ is the solution to the problem

[TABLE]

and $u_{0}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ is the solution to the homogenized problem (2.8).

Proof.

We let $w_{h}:=u_{0}-u_{0}^{h}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ and note that $w_{h}$ is the unique solution to the boundary-value problem

[TABLE]

We recall that $A^{0}\in\mathbb{R}^{n\times n}$ is an elliptic constant matrix. For $h>0$ sufficiently small, by an $H^{2}$ a priori estimate, the Cauchy–Schwarz inequality and Lemma 3.1,

[TABLE]

Finally, we show that for $h>0$ sufficiently small, we have

[TABLE]

with the constant being independent of $h$. This can be seen by rewriting (3.8) as

[TABLE]

Then, again by an $H^{2}$ a priori estimate and Lemma 3.1,

[TABLE]

with constants independent of $h$, i.e., for $h>0$ sufficiently small, (3.9) holds with the constant being independent of $h$. ∎

Finally, we can use an $H^{1}_{0}(\Omega)$-conforming finite element approximation $u_{0}^{h,k}$ to the solution $u_{0}^{h}$ of (3.8), satisfying the error bound

[TABLE]

with constants independent of $h$. By the triangle inequality and the results obtained in this section, we have the following approximation result for $u_{0}$.

Theorem 3.2** ($H^{1}$ approximation of $u_{0}$).**

Assume either that $(\Omega,A,f)\in\mathcal{G}^{0,2}$, or that $(\Omega,A,f)\in\mathcal{H}^{0}$. Then, the approximation $u_{0}^{h,k}$ obtained by the procedure described above satisfies

[TABLE]

Let us now assume either that $(\Omega,A,f)\in\mathcal{G}^{1,2}$ or that $(\Omega,A,f)\in\mathcal{H}^{1}$. Further, assume that for $h>0$ sufficiently small, we have that $u_{0}^{h}\in H^{3}(\Omega)$ with

[TABLE]

where the constant is independent of $h$. The following lemma provides two situations where this is satisfied.

Lemma 3.3**.**

Let $(\Omega,A,f)$ be such that

• $(i)$

$(\Omega,A,f)\in\mathcal{G}^{1,2}$* with $\partial\Omega\in C^{3}$, or*

• $(ii)$

$(\Omega,A,f)\in\mathcal{H}^{1}$* with $\Omega\subset\mathbb{R}^{2}$ being a polygon and $f\in H^{1}_{0}(\Omega)$.*

Then, for $h>0$ sufficiently small, (3.11) holds.

Before we prove Lemma 3.3, we need the following result on the regularity of solutions to Poisson’s problem on convex polygons, see also [13, 15, 16, 17].

Lemma 3.4**.**

Let $\Omega\subset\mathbb{R}^{2}$ be a convex polygonal domain and $f\in H^{1}_{0}(\Omega)$. Then the solution $u\in H^{1}_{0}(\Omega)$ to the problem

[TABLE]

satisfies the bound

[TABLE]

Proof.

First, note that since $\Omega\subset\mathbb{R}^{2}$ is a convex polygonal domain, we have $u\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ with $\|u\|_{H^{2}(\Omega)}\lesssim\|f\|_{L^{2}(\Omega)}$, see [13]. Since $f\in H^{1}_{0}(\Omega)$, there exists a sequence of smooth functions with compact support $(f_{m})_{m}\subset C_{c}^{\infty}(\Omega)$ such that $f_{m}\rightarrow f$ in $H^{1}(\Omega)$. Let $(u_{m})_{m}\subset H^{1}_{0}(\Omega)$ be the sequence of solutions in $H^{1}_{0}(\Omega)$ to $\Delta u_{m}=f_{m}$ in $\Omega$, and note that $(u_{m})_{m}\subset C^{\infty}(\bar{\Omega})$ since the functions $f_{m}$ satisfy compatibility conditions of any order, see [13, Sec. 5.1]. Again we use the $H^{2}$-regularity result for solutions of Poisson’s problem on convex polygons to obtain

[TABLE]

i.e., $u_{m}\rightarrow u$ in $H^{2}(\Omega)$.

Next, we shall use the fact that

[TABLE]

see [17]. We apply (3.13) to the difference of two elements of the sequence $(u_{m})_{m}$ to find that $(u_{m})_{m}$ is a Cauchy sequence in $H^{3}(\Omega)$, using that $f_{m}\rightarrow f$ in $H^{1}(\Omega)$. Thus, $u_{m}\rightarrow u$ in $H^{3}(\Omega)$ and passing to the limit in (3.13) applied to the functions $u_{m}$ yields

[TABLE]

Since $\|u\|_{H^{2}(\Omega)}\lesssim\|f\|_{L^{2}(\Omega)}$, we conclude the bound (3.12). ∎

Remark 3.2**.**

The assumption $f\in H^{1}_{0}(\Omega)$ in Lemma 3.4 can be weakened provided $f$ satisfies certain compatibility conditions, see [13, Theorem 5.1.2.4].

Now we are in a position to prove Lemma 3.3, using standard elliptic regularity theory, Lemma 3.4, and a scaling argument.

Proof of Lemma 3.3.

We start with the case $(i)$. To this end, let $(\Omega,A,f)\in\mathcal{G}^{1,2}$ with $\partial\Omega\in C^{3}$. Then, by elliptic regularity theory, we have $u_{0}^{h}\in H^{3}(\Omega)$. Using elliptic regularity for problem (3.10) yields

[TABLE]

with constants independent of $h$, i.e., for $h>0$ sufficiently small, (3.11) holds with the constant being independent of $h$.

Let us now show the claim for the case $(ii)$. To this end, let $(\Omega,A,f)\in\mathcal{H}^{1}$ with $\Omega\subset\mathbb{R}^{2}$ being a polygon and $f\in H^{1}_{0}(\Omega)$. Since

[TABLE]

is symmetric and elliptic for $h>0$ sufficiently small, there exists an orthogonal matrix $Q_{h}\in\mathbb{R}^{2\times 2}$ with $Q_{h}Q_{h}^{\mathrm{T}}=Q_{h}^{\mathrm{T}}Q_{h}=I_{2}$ such that

[TABLE]

where $\lambda_{h}^{\pm}>0$ are given by

[TABLE]

We note that, by Lemma 3.1, the entries of $B_{h}=(b_{ij}^{h})_{1\leq i,j\leq 2}$ satisfy $b_{ij}^{h}\lesssim h$, and therefore, for $h>0$ sufficiently small, we have $0<\lambda_{h}^{\pm}+(\lambda_{h}^{\pm})^{-1}\lesssim 1$.

The problem (3.8) in the new coordinates reads

[TABLE]

where $U_{h}:=u_{0}^{h}\left(Q_{h}^{\mathrm{T}}\Lambda_{h}^{\frac{1}{2}}\;\cdot\;\right)$, $F_{h}:=f\left(Q_{h}^{\mathrm{T}}\Lambda_{h}^{\frac{1}{2}}\;\cdot\;\right)$, and $P_{h}:=\Lambda_{h}^{-\frac{1}{2}}Q_{h}\Omega$. Note that $P_{h}$ is still a bounded convex polygonal domain and that $F_{h}\in H^{1}_{0}(P_{h})$. By the change of variables formula and the orthogonality of $Q_{h}$,

[TABLE]

Using Lemma 3.4, we have that, for $h>0$ sufficiently small, the solution to (3.14) satisfies

[TABLE]

with constants independent of $h$. It remains to show the bound

[TABLE]

By the change of variables formula and the orthogonality of $Q_{h}$, we obtain similarly as before,

[TABLE]

i.e., we have established the bound (3.15). We conclude that, for $h>0$ sufficiently small, we have (3.11), i.e.,

[TABLE]

where the constant is independent of $h$. ∎

Then an $H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$-conforming finite element approximation $u_{0}^{h,k}$ to the solution $u_{0}^{h}$ of (3.8), that satisfies the error bound

[TABLE]

provides by Lemma 3.2 and the triangle inequality an approximation to $u_{0}$.

Theorem 3.3** ($H^{2}$-norm approximation of $u_{0}$).**

Assume either that $(\Omega,A,f)\in\mathcal{G}^{1,2}$ or that $(\Omega,A,f)\in\mathcal{H}^{1}$, and assume (3.11). Then, the approximation $u_{0}^{h,k}$ obtained by the procedure described above satisfies

[TABLE]

Remark 3.3** (Improvements).**

We note that if we assume that $A\in W^{2,\infty}(Y)$, then we have the following improved results.

• $(i)$

Approximation of $m$: In this case, $m\in H^{2}(Y)$ and we have that

[TABLE]

by choosing $\tilde{v}_{h}=\mathcal{I}_{h}m-\int_{Y}\mathcal{I}_{h}m$, and using an interpolation error bound. Therefore, Theorem 3.1 yields

[TABLE]

• $(ii)$

Approximation of $A^{0}$: By an interpolation error bound and the fact that $m_{h}$ is piecewise linear, one has

[TABLE]

Therefore, the proof of Lemma 3.1 yields

[TABLE]

• $(iii)$

Approximation of $u_{0}$: It follows that the results of Lemma 3.2, Theorem 3.2 and Theorem 3.3 can be improved to second-order convergence in $h$, i.e.,

[TABLE]

for $s=1,2$, respectively.

For the approximation of derivatives of $u_{0}$ of higher than second order, the post-processing method of Babuška in [4] can be used to obtain error bounds in norms involving derivatives of higher order than the energy norm (the norm natural to the problem).

For bounded convex polygonal domains $\Omega\subset\mathbb{R}^{2}$, an $H^{2}$-conforming approximation to the solution of (3.8) can be obtained as follows. Assume that $f\in H^{1}_{0}(\Omega)$ so that (3.11) holds. Consider a shape-regular triangulation of $\Omega$ into triangles with longest edge $k>0$, and let

[TABLE]

be an appropriate finite element space. In practice, the Hsieh–Clough–Tocher element and the Argyris element can be used as $H^{2}$-conforming elements. Then, for $h>0$ sufficiently small, standard finite element analysis can be used to show that there is a unique function $u_{0}^{h,k}\in V_{k}$ such that

[TABLE]

and that the error bound (3.16) holds.

3.2. Approximation of the Corrector

We now address problem (2.9) and present a method for $A\in W^{2,\infty}(Y)$. To simplify the notation and the arguments, we assume that we know the invariant measure $m$ and the matrix $A^{0}=(a_{ij}^{0})_{1\leq i,j\leq n}$ exactly instead of working with our approximation $A_{h}^{0}$.

For a given $Y$-periodic right-hand side $g\in W^{2,\infty}(Y)$, we address the problem

[TABLE]

Obtaining an approximation for second-order derivatives via finite elements is not straightforward since the natural solution space is $W_{\text{per}}(Y)$. We present a method of successively approximating higher derivatives.

Let $\chi_{h}$ be a $W_{\text{per}}(Y)$-conforming finite element approximation to $\chi$, i.e.,

[TABLE]

with $V_{h}\subset W_{\text{per}}(Y)$ finite-dimensional, and satisfying the error estimate

[TABLE]

Let $r\in\{1,\dots,n\}$ and write $\xi_{r}:=\partial_{r}\chi$. Then, using the equation

[TABLE]

we find that weakly, there holds

[TABLE]

Further, we claim that $\xi_{r}\in W_{\text{per}}(Y)$. Indeed, this follows from the regularity and periodicity of $\chi$ and

[TABLE]

Therefore, $\xi_{r}\in W_{\text{per}}(Y)$ satisfies

[TABLE]

Now we use our $H^{1}$-conforming approximation for $\chi$ for the right-hand side and use a $W_{\text{per}}(Y)$-conforming finite element method for approximating the solution $v\in W_{\text{per}}(Y)$ to the following problem:

[TABLE]

where $c$ is such that this problem admits a unique solution (such that the solvability condition (2.5) is satisfied). By looking at the problem for $v-\xi_{r}$, one obtains the comparison result

[TABLE]

Let $v_{h}$ be a $W_{\text{per}}(Y)$-conforming finite element approximation of (3.18) satisfying

[TABLE]

for some constant $C=C(\|A\|_{W^{2,\infty}(Y)})>0$. Then, using the triangle inequality, we obtain

[TABLE]

for some constant $C=C(\|A\|_{W^{2,\infty}(Y)})>0$. Using this procedure for $r=1,\dots,n$, we eventually obtain approximations to derivatives of order up to two of $\chi$.

3.3. Approximation of $u_{\varepsilon}$

We assume either that $(\Omega,A,f)\in\mathcal{G}^{2,2}$ or that $(\Omega,A,f)\in\mathcal{H}^{2}$. Let $n\in\{2,3\}$, $\varepsilon\in(0,1]$, and assume that

[TABLE]

Then we know that (2.11) is satisfied, and by Theorem 2.5 we have that

[TABLE]

where $u_{0}$ is the solution to the homogenized problem, and $\chi_{ij}$ are the corrector functions given as the solutions to (2.9). This result can be used to construct an approximation of $u_{\varepsilon}$, i.e., to the solution of problem (2.2) for small $\varepsilon$. We note that (3.19) implies that

[TABLE]

This leads to the following approximation result for $u_{\varepsilon}$.

Theorem 3.4** (Approximation of $u_{\varepsilon}$).**

In the situation described above, let $(u_{0,h})_{h>0}\subset H^{2}(\Omega)$ be a family of $H^{2}$-conforming approximations for $u_{0}$ satisfying the error bound

[TABLE]

and for $1\leq i,j,k,l\leq n$, let $(z_{ij,h}^{kl})_{h>0}\subset L^{2}_{\text{per}}(Y)$ be a family of $L^{2}$ approximations for $\partial_{kl}^{2}\chi_{ij}$ satisfying the error bound

[TABLE]

Then, by writing

[TABLE]

we have that

[TABLE]

Proof.

We use (3.20) and the triangle inequality to obtain

[TABLE]

and for $1\leq k,l\leq n$,

[TABLE]

It remains to study the last term on the right-hand side of the above inequality. For fixed $1\leq i,j\leq n$, we use again the triangle inequality to obtain

[TABLE]

In the last step, we used that by the transformation formula and periodicity (cover $\Omega/\varepsilon$ by $\mathcal{O}(\varepsilon^{-n})$ many cells of unit length), there holds

[TABLE]

We claim that

[TABLE]

Indeed, we use the triangle inequality, (3.21) and the fact that $\chi_{ij}\in W^{2,\infty}(Y)$ to obtain

[TABLE]

∎

The approximations of $u_{0}$ and the corrector functions can be obtained as described in Section 3.1 and 3.2. Let us conclude this section by remarking that if the second derivatives of the corrector functions are approximated in the space $L^{\infty}(Y)$ or if the solution to the homogenized problem is approximated in the space $W^{2,\infty}(\Omega)$, then one obtains by a similar proof an approximation result for the second derivatives of $u_{\varepsilon}$ in $L^{2}(\Omega)$.

Remark 3.4**.**

If $(z_{ij,h}^{kl})_{h>0}\subset L^{\infty}_{\text{per}}(Y)$ is a family of $L^{\infty}$ approximations for $\partial_{kl}^{2}\chi_{ij}$ satisfying the error bound

[TABLE]

and $(u_{0,h})_{h>0}$ is as in Theorem 3.4, then we have that

[TABLE]

The same holds true when $(u_{0,h})_{h>0}\subset W^{2,\infty}(\Omega)$ is a family of $W^{2,\infty}$-conforming approximations for $u_{0}$ satisfying the error bound

[TABLE]

and $(z_{ij,h}^{kl})_{h>0}$ is as in Theorem 3.4.

3.4. Nonuniformly Oscillating Coefficients

In this section, we discuss the case of nonuniformly oscillating coefficients, i.e., coefficients depending on $x$ and $\frac{x}{\varepsilon}$. We consider the problem

[TABLE]

where the triple $(\Omega,A,f)$ satisfies one of the following sets of assumptions.

Definition 3.1** (Sets of assumptions $\mathcal{G},\mathcal{H}$).**

For $m\in\mathbb{N}_{0}$, we write

• $(i)$

$(\Omega,A,f)\in\mathcal{G}$* if and only if $\Omega\subset\mathbb{R}^{n}$ is a bounded $C^{2,\gamma}$ domain, $f\in L^{2}(\Omega)$, and $A=A^{\mathrm{T}}:\Omega\times\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ satisfies*

[TABLE]

• $(ii)$

$(\Omega,A,f)\in\mathcal{H}$* if and only if $\Omega\subset\mathbb{R}^{n}$ is a bounded convex domain, $f\in L^{2}(\Omega)$, and $A=A^{\mathrm{T}}:\Omega\times\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ satisfies (3.23) and*

[TABLE]

In view of Remark 2.1, we see that the Cordes condition (3.24) is always satisfied when $n=2$. Well-posedness to the problem (3.22) is guaranteed by the following theorem, see [12, Theorem 9.15] and [19, Theorem 3].

Theorem 3.5** (Existence and uniqueness of strong solutions).**

Assume either that $(\Omega,A,f)\in\mathcal{G}$, or that $(\Omega,A,f)\in\mathcal{H}$. Then, for any $\varepsilon>0$, the problem (3.22) admits a unique solution $u_{\varepsilon}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$.

As in Section 2, uniform a priori estimates for the solution to (3.22) allow passage to the limit in equation (3.22), see [5, 6]. The coefficient matrix of the homogenized problem now depends on the slow variable $x$, and is obtained by integrating against an invariant measure. Corrector results can then be shown as before.

Theorem 3.6** (Nonuniformly oscillating coefficients).**

Assume that $\varepsilon\in(0,1]$ and either that $(\Omega,A,f)\in\mathcal{G}$, or that $(\Omega,A,f)\in\mathcal{H}$. Then the following assertions hold.

• $(i)$

Uniform a priori estimate: The solution $u_{\varepsilon}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ to (3.22) satisfies

[TABLE]

• $(ii)$

Homogenization: The solution $u_{\varepsilon}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ to (3.22) converges weakly in $H^{2}(\Omega)$ to the solution $u_{0}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ of the homogenized problem

[TABLE]

with $A^{0}:\Omega\rightarrow\mathbb{R}^{n\times n}$ given by

[TABLE]

where $m=m(x,y)$ is the unique function $m:\Omega\times\mathbb{R}^{n}\rightarrow\mathbb{R}$ with $m\in C(\bar{\Omega}\times\mathbb{R}^{n})$, $0<\bar{m}\leq m\leq M$ for some constants $\bar{m},M>0$, such that

[TABLE]

for any fixed $x\in\Omega$. The function $m$ is called the invariant measure.

• $(iii)$

Corrector estimate: Assume that $f\in H^{2}(\Omega)$ and $u_{0}\in H^{4}(\Omega)\cap W^{2,\infty}(\Omega)$. Introducing the corrector function $\chi_{ij}$, $1\leq i,j\leq n$, as the solution to

[TABLE]

we have that

[TABLE]

Proof.

$(i)$ For $(\Omega,A,f)\in\mathcal{H}$, one shows similarly to the proof of [19, Theorem 3] and Theorem 2.2 that

[TABLE]

For $(\Omega,A,f)\in\mathcal{G}$, the claim follows from the method of freezing coefficients, using the uniform estimate from Theorem 2.2 for the operators $L_{x_{0}}:=A\left(x_{0},\frac{\cdot}{\varepsilon}\right):D^{2}$ for fixed $x_{0}\in\Omega$.

$(ii)$ The uniform estimate from $(i)$ yields weak convergence in $H^{2}(\Omega)$ and strong convergence in $H^{1}(\Omega)$ for a subsequence of $(u_{\varepsilon})_{\varepsilon>0}$ to some limit function $u_{0}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$. We multiply (3.22) by $m\left(\cdot,\frac{\cdot}{\varepsilon}\right)$ and follow the transformation performed in [5] to find that

[TABLE]

holds weakly, where $\tilde{A}:=Am$ and $v^{\varepsilon}$ denotes $v\left(\cdot,\frac{\cdot}{\varepsilon}\right)$. Passing to the limit, we obtain that $u_{0}$ is a weak solution of (3.25). We conclude the proof by noting that (3.25) admits a unique strong solution, since $A^{0}$ is uniformly elliptic and Lipschitz continuous on $\bar{\Omega}$, see [12, 13].

$(iii)$ This can be proved similarly to Theorem 2.4 and Theorem 2.5, using that, by the assumptions made on $A$ and elliptic regularity, we have

[TABLE]

for any $1\leq i,j,k,l\leq n$. ∎

Let us explain how the numerical scheme from Section 3.1 can be used for the numerical homogenization of (3.22).

First, we consider a triangulation $\mathcal{T}_{k}$ on $\bar{\Omega}$ consisting of nodes $\{x_{i}\}_{i\in I}$ with grid size $k>0$, and a triangulation $\mathcal{T}_{h}$ on $Y$ with grid size $h>0$. Then, for any $i\in I$, we can use the scheme from Section 3.1 (see Theorem 3.1) to obtain an approximation $m^{i}_{h}\in H^{1}(Y)$ to $m_{x_{i}}=m(x_{i},\cdot)$ such that

[TABLE]

Further, we obtain that

[TABLE]

is an approximation to $A^{0}(x_{i})$ (see Lemma 3.1),

[TABLE]

Now we define $A_{h,k}^{0}$ to be a continuous piecewise linear function on $\Omega$ such that

[TABLE]

Then, using (3.26) and denoting the continuous piecewise linear interpolant of a function $\phi$ on the grid $\{x_{i}\}_{i\in I}$ by $\mathcal{I}_{k}\phi$, we have

[TABLE]

We observe that, similarly to the proof of Lemma 3.2, we obtain that the solution $u_{0}^{h,k}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ to

[TABLE]

satisfies, for $h,k>0$ sufficiently small,

[TABLE]

and in view of (3.27),

[TABLE]

where $u_{0}$ is the solution to the homogenized problem (3.25). Finally, the solution to (3.28) can be approximated by a standard finite element method on the triangulation $\mathcal{T}_{k}$ which yields an approximation $u_{0,h,k}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ to $u_{0}$ in the $H^{2}(\Omega)$-norm.

The approximation of $u_{\varepsilon}$ can be obtained based on the corrector estimate from Theorem 3.6 analogously as in Section 3.3.

4. Numerical Experiments

4.1. Problem with a Known $u_{0}$

We consider the homogenization problem

[TABLE]

on the domain

[TABLE]

with the matrix-valued map

[TABLE]

and the right-hand side $f:\Omega\rightarrow\mathbb{R}$ to be specified below. We observe that the matrix-valued function $A$ satisfies (2.1) with $q=\infty$. Further, note that

[TABLE]

depends only on the first coordinate of $y=(y_{1},y_{2})\in\mathbb{R}^{2}$; see Figure 1.

In this case we know that the homogenized problem is given by

[TABLE]

where $A^{0}\in\mathbb{R}^{2\times 2}$ denotes the constant matrix

[TABLE]

with $m$ being the invariant measure

[TABLE]

see [11]. Explicit computation yields that

[TABLE]

We consider the right-hand side given by

[TABLE]

Then it is straightforward to check that the exact solution $u_{0}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ to the homogenized problem (4.2) is given by

[TABLE]

Note that we are in the situation $(\Omega,A,f)\in\mathcal{H}^{2}$, that $f=0$ in the corners of $\Omega$ and that $u_{0}\in H^{4}(\Omega)$.

We use the scheme presented in Section 3.1 to approximate $m$, $A^{0}$ and $u_{0}$. We use the same mesh for approximating $m$ and $u_{0}$. The Hsieh–Clough–Tocher (HCT) element in FreeFem++ is used in the formulation (3.17) for the $H^{2}$ approximation of $u_{0}$; see [14]. The gradient on the boundary is set to be the gradient of an $H^{1}$ approximation by $\mathbb{P}_{2}$ elements on a fine mesh.

Figure 2 shows the error in the approximation of $m$ and $A^{0}$. For the approximation of the invariant measure we observe convergence of order

[TABLE]

and superconvergence of order $\mathcal{O}(h^{2})$ for $h>0$ when grid points fall on the line $\{y_{1}=\frac{1}{2}\}$, which is the set along which $\partial_{1}m$ possesses a jump. The observed rate of convergence (4.3) is consistent with Theorem 3.1. Indeed, we have $m\in H^{\frac{3}{2}-\tilde{\varepsilon}}(Y)$ for any $\tilde{\varepsilon}>0$, and Theorem 3.1 yields

[TABLE]

by making the choice $\tilde{v}_{h}=\mathcal{I}_{h}m-\int_{Y}\mathcal{I}_{h}m$, and using an interpolation error bound. In connection with the superconvergence we note that $\left.m\right\rvert_{(0,\frac{1}{2})\times(0,1)}\in H^{2}((0,\frac{1}{2})\times(0,1))$ and $\left.m\right\rvert_{(\frac{1}{2},1)\times(0,1)}\in H^{2}((\frac{1}{2},1)\times(0,1))$. For the approximation of the matrix $A^{0}$, we observe second-order convergence.

Concerning the approximation of $u_{\varepsilon}$, from Sections 2 and 3.3 we obtain that

[TABLE]

where $\chi_{ij}$ ($1\leq i,j\leq 2$) denotes the solution to

[TABLE]

Note that since $A$ only depends on $y_{1}$, we have that

[TABLE]

Therefore, there holds

[TABLE]

For the numerical approximation, we replace $u_{\varepsilon}$ by an $H^{2}$-conforming finite element approximation on a fine mesh, based on the formulation

[TABLE]

where $H:=H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$. To this end, we use again the HCT element and set the gradient on the boundary to be the gradient of an $H^{1}$ approximation by $\mathbb{P}_{2}$ elements on a fine mesh.

Figure 3 shows the error in the approximation of $u_{0}$ and we observe second-order convergence. Further, with the exact $u_{0}$ being available, we can compute the error (4.4) for different values of $\varepsilon$; see Figure 3. We observe first-order convergence as $\varepsilon$ tends to zero, as expected from (4.4).

4.2. Problem with an Unknown $u_{0}$

Next, let us consider the problem (4.1) with the same domain $\Omega$ and matrix-valued function $A$ as before, but with the right-hand side given by

[TABLE]

Note that we are in the situation $(\Omega,A,f)\in\mathcal{H}^{2}$. Further, since the right-hand side $f\in H^{2}(\Omega)$ of the homogenized problem (4.2) satisfies $f=0$ at the corners of $\Omega$, the solution $u_{0}$ to (4.2) belongs to the class $H^{4}(\Omega)$; see [16, Prop. 2.6].

As before, we use the scheme presented in Section 3.1 to approximate $m$, $A^{0}$ and $u_{0}$. Using the second-order $H^{2}$ approximation $u_{0,h}$ to $u_{0}$ obtained as previously described,

[TABLE]

we have that

[TABLE]

Figure 4 shows the squared error (4.5) of the approximation of $u_{\varepsilon}$ for different grid sizes and $\varepsilon=\frac{1}{100}$ fixed. We observe fourth-order convergence in $h$ for the squared error as expected from (4.5).

4.3. Nonuniformly Oscillating Coefficients

We consider the homogenization problem

[TABLE]

on the domain

[TABLE]

with the matrix-valued map $A:\Omega\times\mathbb{R}^{2}\rightarrow\mathbb{R}^{2\times 2}$,

[TABLE]

and the right-hand side $f:\Omega\rightarrow\mathbb{R}$ to be specified below. We observe that the matrix-valued function $A$ satisfies (3.23) with $q=\infty$. Further, note that it is of the form

[TABLE]

In this case we know that the homogenized problem is given by

[TABLE]

where $A^{0}:\Omega\rightarrow\mathbb{R}^{2\times 2}$ is given by

[TABLE]

with $m$ being the invariant measure

[TABLE]

see [11]. Therefore, we have

[TABLE]

We consider the right-hand side given by

[TABLE]

Then it is straightforward to check that the exact solution $u_{0}\in H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$ to the homogenized problem (4.7) is given by

[TABLE]

Note that the assumptions of Theorem 3.6 $(iii)$ are satisfied.

For $k>0$ such that $\frac{1}{k}\in\mathbb{N}$, we take a triangulation $\mathcal{T}_{k}$ on $\bar{\Omega}$ consisting of nodes $\{(sk,rk)\}_{s,r=0,\dots,1/k}$, and a triangulation $\mathcal{T}_{h}$ on $Y$ with grid size $h=\frac{k}{4}$. We use the scheme presented in Section 3.4 to approximate $A^{0}$ and $u_{0}$, and we observe second-order convergence; see Figure 5.

For the approximation of $u_{\varepsilon}$, Theorem 3.6 yields

[TABLE]

where $\chi_{ij}$ ($1\leq i,j\leq 2$) denotes the solution to

[TABLE]

We observe that we have

[TABLE]

Therefore, we have that

[TABLE]

where $\left[\partial_{y_{i}y_{i}}^{2}\chi_{ii}\right]^{\varepsilon}:=\left(\partial_{y_{i}y_{i}}^{2}\chi_{ii}\right)\left(\cdot,\frac{\cdot}{\varepsilon}\right)$ for $i\in\{1,2\}$. For the numerical approximation, we replace $u_{\varepsilon}$ by an $H^{2}$-conforming finite element method on a fine mesh, based on the formulation

[TABLE]

where $H:=H^{2}(\Omega)\cap H^{1}_{0}(\Omega)$. To this end, we use again the HCT element and set the gradient on the boundary to be the gradient of an $H^{1}$ approximation by $\mathbb{P}_{2}$ elements on a fine mesh.

Finally, let us consider the problem (4.6) with the same domain $\Omega$ and matrix-valued function $A$ as before, but with the right-hand side given by

[TABLE]

Note that we are in the situation $(\Omega,A,f)\in\mathcal{H}$. Further, since the right-hand side $f\in H^{2}(\Omega)$ of the homogenized problem (4.7) satisfies $f=0$ at the corners of $\Omega$, the solution $u_{0}$ to (4.7) belongs to the class $H^{4}(\Omega)$, see [16, Prop. 2.6].

Using the second-order $H^{2}$-conforming approximation $u_{0,k}$ to $u_{0}$ obtained as previously described (again with $h=\frac{k}{4}$),

[TABLE]

we have that

[TABLE]

Figure 6 shows the squared error (4.9) of the approximation of $u_{\varepsilon}$ for different grid sizes and $\varepsilon=\frac{1}{50}$ fixed. We observe fourth-order convergence in $k$ for the squared error as expected from (4.9).

5. Conclusion

In this paper we introduced a scheme for the numerical approximation of elliptic problems in nondivergence-form with rapidly oscillating coefficients on $C^{2,\gamma}$ and polygonal domains, which is based on a $W^{2,p}$ corrector estimate for such problems derived in the first part of this work.

We proved an optimal-order error bound for a finite element approximation of the corresponding invariant measure using continuous $Y$-periodic piecewise linear basis functions on a shape-regular triangulation of the unit cell $Y$ under weak regularity assumptions on the coefficients. The coefficients are integrated against the so obtained approximation of the invariant measure after piecewise linear interpolation on the mesh to obtain an approximation of the constant coefficient-matrix of the homogenized problem. Using an $H^{2}$ comparison result for the solution of this perturbed problem, we eventually obtained an approximation of the solution $u_{0}$ to the homogenized problem in the $H^{2}$-norm. In the case of a polygonal domain in two space dimensions, we made use of compatibility conditions for the source term to ensure sufficiently high Sobolev-regularity of $u_{0}$.

We obtained an approximation to the solution $u_{\varepsilon}$ of the original problem, i.e., the problem with oscillating coefficients, by making use of the $H^{2}$ approximation of $u_{0}$, finite element approximations to second-order derivatives of the corrector functions, as well as an $H^{2}$ corrector result. A method of successively approximating higher derivatives for the approximation of corrector functions was provided and analyzed. The corrector functions are necessary in order to obtain an approximation of $D^{2}u_{\varepsilon}$ whereas the task of approximating $u_{\varepsilon}$ in the $H^{1}$-norm can be achieved using only an $H^{1}$ approximation of $u_{0}$.

Furthermore, we generalized our results to the case of nonuniformly oscillating coefficients, i.e., we derived an analogous corrector result and studied the approximation of the solution $u_{0}$ to the homogenized problem and the solution $u_{\varepsilon}$ of the $\varepsilon$-dependent problem in this case.

In the final part of the paper, we presented numerical experiments matching the theoretical results for problems with both known and unknown $u_{0}$, as well as problems with nonuniformly oscillating coefficients. We illustrated the performance of the scheme for the approximation of the invariant measure, the solution $u_{0}$ to the homogenized problem and the solution $u_{\varepsilon}$ to the problem involving oscillating coefficients for a fixed value of $\varepsilon$.

Future work will focus on weakening of the regularity assumptions on the coefficients and the approximation of fully nonlinear nondivergence-form problems with oscillating coefficients such as the Hamilton–Jacobi–Bellman equation.

Acknowledgements

This work was supported by the UK Engineering and Physical Sciences Research Council [EP/L015811/1].