Strong convergence of a linearization method for semi-linear elliptic equations with variable scaled production

TL;DR

This paper introduces a linearization algorithm for semi-linear elliptic equations with variable scaled production in perforated media, ensuring strong convergence and solvability at micro and macro scales, validated through finite element simulations.

Contribution

The work develops a novel linearization scheme that guarantees strong $H^1$ convergence for complex elliptic problems with degenerate production, combining micro-scale contraction and macro-scale homogenization.

Findings

Proves $H^1$-contraction at the micro-scale using energy methods.

Establishes convergence at the macro-scale via homogenization and corrector estimates.

Numerical results confirm the efficiency and convergence of the proposed algorithm.

Abstract

This work is devoted to the development and analysis of a linearization algorithm for microscopic elliptic equations, with scaled degenerate production, posed in a perforated medium and constrained by the homogeneous Neumann-Dirichlet boundary conditions. This technique plays two roles: to guarantee the unique weak solvability of the microscopic problem and to provide a fine approximation in the macroscopic setting. The scheme systematically relies on the choice of a stabilization parameter in such a way as to guarantee the strong convergence in $H^1$ norm for both the microscopic and macroscopic problems. In the standard variational setting, we prove the $H^1$-type contraction at the micro-scale based on the energy method. Meanwhile, we adopt the classical homogenization result in line with corrector estimate to show the convergence of the scheme at the macro-scale. In the numerical section, we use the standard finite element method to assess the efficiency and convergence of our proposed algorithm.

Tables

Table 1. Table 1 : The errors: E 1 subscript 𝐸 1 E_{1} between u ε subscript 𝑢 𝜀 u_{\varepsilon} and u 0 k superscript subscript 𝑢 0 𝑘 u_{0}^{k} , and between ∇ u ε ∇ subscript 𝑢 𝜀 \nabla u_{\varepsilon} and ∇ u 0 k ∇ superscript subscript 𝑢 0 𝑘 \nabla u_{0}^{k} ; the relative error E ​ 2 𝐸 2 E2 between u ε subscript 𝑢 𝜀 u_{\varepsilon} and u 0 k superscript subscript 𝑢 0 𝑘 u_{0}^{k} (with k = 2 ) k=2)

$\epsilon$	0.5	0.25	$0.166$	0.1	0.083	0.05	0.025
$E_1(u_{\epsilon ,h},u_{0,h}^k)$	0.0113	0.0040	0.0026	0.0018	0.0017	0.0016	0.00158
$E_1(\nabla u_{\epsilon ,h},\nabla u_{0,h}^k)$	0.2167	0.1857	0.1824	0.1761	0.1755	0.1746	0.17425
$E_2(u_{\epsilon ,h},u_{0,h}^k)$	0.1854	0.0562	0.0357	0.0242	0.0230	0.0217	0.02128
${\max }_h$	0.04	0.02	0.014	0.008	0.007	0.006	0.0058

Table 2. Table 2 : The relative errors E 2 subscript 𝐸 2 E_{2} between u ε , h k superscript subscript 𝑢 𝜀 ℎ 𝑘 u_{\varepsilon,h}^{k} and u 0 , h k subscript 𝑢 0 superscript ℎ 𝑘 u_{0,h^{k}} , between u ε , h subscript 𝑢 𝜀 ℎ u_{\varepsilon,h} and u ε , h k superscript subscript 𝑢 𝜀 ℎ 𝑘 u_{\varepsilon,h}^{k} , and between ∇ u ε , h ∇ subscript 𝑢 𝜀 ℎ \nabla u_{\varepsilon,h} and ∇ u ε , h k ∇ superscript subscript 𝑢 𝜀 ℎ 𝑘 \nabla u_{\varepsilon,h}^{k} for the fixed value of ε = 0.25 𝜀 0.25 \varepsilon=0.25 as k 𝑘 k increases up to 4

$k$	1	2	3	4
$E_2(u_{\epsilon ,h}^k,u_{0,h}^k)$	0.06591	0.06441	0.06413	0.06408
$E_2(u_{\epsilon ,h},u_{\epsilon ,h}^k)$	0.13964	0.01764	0.00219	0.00027
$E_2(\nabla u_{\epsilon ,h},\nabla u_{\epsilon ,h}^k)$	0.13623	0.01720	0.00214	0.00027

Equations

$$(P_{\epsilon }):\{\begin{array}{l}{L}_{\epsilon }{u}_{\epsilon }+{\epsilon }^{\alpha }\mathcal{R}(u_{\epsilon })=f(x)\text{in}{\Omega }^{\epsilon },\alpha \ge 0,\\ -\text{A}(x/\epsilon )\nabla u_{\epsilon }\cdot n=0\text{across}{\Gamma }^{\epsilon },\\ u_{\epsilon }=0\text{across}{\Gamma }^{\text{ext}},\end{array}$$

$$L_{\epsilon }u=\nabla \cdot \left(-\text{A}\left(\frac{x}{\epsilon }\right)\nabla u\right)=\sum _{i,j=1}^d\frac{\partial }{\partial x_i}\left(-a_{ij}\left(\frac{x}{\epsilon }\right)\frac{\partial u}{\partial x_j}\right).$$

$$\displaystyle\begin{array}[]{ccccccc}\left(P_{\varepsilon}\right)&\longrightarrow&\left(P_{\varepsilon}^{k}\right)&\longrightarrow&\left(P_{0}^{k}\right)\\ \downarrow&&\downarrow&&\downarrow&\\ u_{\varepsilon}&\longrightarrow&u_{\varepsilon}^{k}&\longrightarrow&u_{0}^{k}\end{array}$$

$$\displaystyle\begin{array}[]{ccccccc}\left(P_{\varepsilon}\right)&\longrightarrow&\left(P_{0}\right)&\longrightarrow&\left(P_{0}^{k}\right)\\ \downarrow&&\downarrow&&\downarrow&\\ u_{\varepsilon}&\longrightarrow&u_{0}&\longrightarrow&u_{0}^{k}\end{array}$$

$$0<{\delta }_0\le F^{\prime }\left(x\right)\le {\delta }_1\text{for a.e. }x\in \mathbb{R}\backslash \mathcal{B}\left(x_0,r_{x_0}\right),$$

$$V^{\epsilon }:=\left\{u\in H^1\left({\Omega }^{\epsilon }\right):u=0\text{on }{\Gamma }^{\text{ext}}\right\},$$

$${\Vert u\Vert }_{V^{\epsilon }}={\left({\int }_{{\Omega }^{\epsilon }}{\left|\nabla u\right|}^{2}dx\right)}^{1/2}.$$

$${\Vert u\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2\le C_p{\Vert \nabla u\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2\text{for any }u\in V^{\epsilon }.$$

$$\underline{\gamma }{\left|\xi \right|}^2\le \sum _{i,j=1}^d{a}_{ij}\left(\frac{x}{\epsilon }\right){\xi }_i{\xi }_j\le \overline{\gamma }{\left|\xi \right|}^2\text{for all }\xi \in {\mathbb{R}}^d\text{and }1\le i,j\le d.$$

$$a(u_{\epsilon },\phi )+{\epsilon }^{\alpha }⟨\mathcal{R}(u_{\epsilon }),\phi {⟩}_{L^2({\Omega }^{\epsilon })}=⟨f,\phi {⟩}_{L^2({\Omega }^{\epsilon })},$$

$$a(u,\phi ):={\int }_{{\Omega }^{\epsilon }}\text{A}(x/\epsilon )\nabla u\cdot \nabla \phi dx.$$

$$\underset{x\in \left[-1,1\right]}{\sup }\left|\bar{f}-p_n\right|\le C{n}^{-k}\omega \left({\bar{f}}^{(k)},n^{-1}\right),$$

$$a(u_{\epsilon }^k,\phi )+M⟨u_{\epsilon }^k,\phi {⟩}_{L^2({\Omega }^{\epsilon })}=⟨f,\phi {⟩}_{L^2({\Omega }^{\epsilon })}+M⟨u_{\epsilon }^{k-1},\phi {⟩}_{L^2({\Omega }^{\epsilon })}-{\epsilon }^{\alpha }⟨{\mathcal{R}}_{{\gamma }_k}(u_{\epsilon }^{k-1}),\phi {⟩}_{L^2({\Omega }^{\epsilon })},$$

$$a\left(u_{\epsilon }^1,\phi \right)+M{⟨u_{\epsilon }^1,\phi ⟩}_{L^2\left({\Omega }^{\epsilon }\right)}={⟨\tilde{f},\phi ⟩}_{L^2\left({\Omega }^{\epsilon }\right)},$$

$$p_k+q_k\le a_k+b_k{q}_{k-1},k\ge 2,$$

$$p_k+q_k\le a_k+\sum _{j=2}^{k-1}{a}_j\prod _{i=j+1}^k{b}_i+q_1\prod _{i=2}^k{b}_i,k\ge 3.$$

$$\frac{\underline{\gamma }}{\eta +{\delta }_1+\underline{\gamma }{C}_p^{-1}}{\Vert \nabla u_{\epsilon }^k-\nabla u_{\epsilon }^{k-1}\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2+{\Vert u_{\epsilon }^k-u_{\epsilon }^{k-1}\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2$$

$$\le \frac{C{\gamma }_{k-1}^{2\sigma }}{{\gamma }_k}+C\sum _{j=2}^{k-1}{\left(\frac{M}{\underline{\gamma }{C}_p^{-1}+M}\right)}^{k-j}\frac{{\gamma }_{j-1}^{2\sigma }}{{\gamma }_j}+{\Vert u_{\epsilon }^1\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2{\left(\frac{M}{\underline{\gamma }{C}_p^{-1}+M}\right)}^{k-1}.$$

$$a(w_{\epsilon }^k,\phi )+M⟨w_{\epsilon }^k,\phi {⟩}_{L^2({\Omega }^{\epsilon })}=M⟨w_{\epsilon }^{k-1},\phi {⟩}_{L^2({\Omega }^{\epsilon })}-{\epsilon }^{\alpha }⟨{\mathcal{R}}_{{\gamma }_k}(u_{\epsilon }^{k-1})-{\mathcal{R}}_{{\gamma }_{k-1}}(u_{\epsilon }^{k-2}),\phi {⟩}_{L^2({\Omega }^{\epsilon })}.$$

$${\int }_{{\Omega }^{\epsilon }}\mathbf{A}\left(\frac{x}{\epsilon }\right){\left|\nabla w_{\epsilon }^k\right|}^{2}dx+M{\Vert w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2$$

$$={\epsilon }^{\alpha }{⟨{\mathcal{R}}_{{\gamma }_{k-1}}\left(u_{\epsilon }^{k-2}\right)-{\mathcal{R}}_{{\gamma }_k}\left(u_{\epsilon }^{k-2}\right),w_{\epsilon }^k⟩}_{L^2\left({\Omega }^{\epsilon }\right)}+{⟨g_{{\gamma }_k}\left(u_{\epsilon }^{k-2}\right)-g_{{\gamma }_k}\left(u_{\epsilon }^{k-1}\right),w_{\epsilon }^k⟩}_{L^2\left({\Omega }^{\epsilon }\right)}.$$

$$\left|{\mathcal{R}}_{{\gamma }_{k-1}}-{\mathcal{R}}_{{\gamma }_k}\right|\le C\left({\gamma }_{k-1}^{\sigma }+{\gamma }_k^{\sigma }\right),$$

$$\underline{\gamma }{\Vert \nabla w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2+M{\Vert w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2$$

$$\le \underset{:=I_1}{\underbrace{C{\epsilon }^{\alpha }\left({\gamma }_{k-1}^{\sigma }+{\gamma }_k^{\sigma }\right){\Vert w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}}}+\underset{:=I_2}{\underbrace{\left(M-{\epsilon }^{\alpha }{\gamma }_k{\delta }_0\right){\Vert w_{\epsilon }^{k-1}\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}{\Vert w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}}}.$$

$$I_1\le \frac{C{\epsilon }^{\alpha }{\gamma }_{k-1}^{2\sigma }}{{\delta }_0{\gamma }_k}+\frac{{\epsilon }^{\alpha }{\delta }_0{\gamma }_k}{2}{\Vert w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2.$$

$$I_2\le \frac{M-{\epsilon }^{\alpha }{\gamma }_k{\delta }_0}{2}\left({\Vert w_{\epsilon }^{k-1}\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2+{\Vert w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2\right).$$

$$2\underline{\gamma }{\Vert \nabla w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2+M{\Vert w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2\le \frac{C{\gamma }_{k-1}^{2\sigma }}{{\delta }_0{\gamma }_k}{\epsilon }^{\alpha }+\left(M-{\epsilon }^{\alpha }{\gamma }_k{\delta }_0\right){\Vert w_{\epsilon }^{k-1}\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2.$$

$$\frac{\underline{\gamma }}{\eta +{\delta }_1+\underline{\gamma }{C}_p^{-1}}{\Vert \nabla w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2+{\Vert w_{\epsilon }^k\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2\le \frac{C{\gamma }_{k-1}^{2\sigma }}{{\gamma }_k}{\epsilon }^{\alpha }+\frac{M}{\underline{\gamma }{C}_p^{-1}+M}{\Vert w_{\epsilon }^{k-1}\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2,$$

$$\frac{\underline{\gamma }}{\eta +{\delta }_1+\underline{\gamma }{C}_p^{-1}}$$

$$\le a_k+\sum _{j=2}^{k-1}{a}_j\prod _{i=j+1}^k{b}_i+{\Vert u_{\epsilon }^1\Vert }_{L^2\left({\Omega }^{\epsilon }\right)}^2\prod _{i=2}^k{b}_i,k\ge 3,$$

$$a_k:=\frac{C{\gamma }_{k-1}^{2\sigma }}{{\gamma }_k},b_k:=\frac{M}{\underline{\gamma }{C}_p^{-1}+M}.$$

$$\prod _{i=j+1}^k{b}_i={\left(\frac{M}{\underline{\gamma }{C}_p^{-1}+M}\right)}^{k-j}\text{for }j\ge 1.$$

Full text

1 Introduction

1.1 Microscopic problem

Let $\Omega^{\varepsilon}$ be a Lipschitz perforated domain contained in a polygonal bounded domain $\Omega\subset\mathbb{R}^{d}$ $(d=2,3)$. In this sense, $\Omega^{\varepsilon}$ possesses a uniformly periodic microstructure defined by a length scale $\varepsilon$. This $\varepsilon$ is a small parameter $(0<\varepsilon\ll 1)$ since the size of the pores are usually much smaller than the characteristic length of the reservoir. We are herein concerned with the asymptotic behavior in a stationary case of the function $u_{\varepsilon}:\Omega^{\varepsilon}\to\mathbb{R}$ that describes the spread of concentration of solutes dissolved in a saturated porous tissue shaped by the perforated domain $\Omega^{\varepsilon}$ with a cubic periodicity cell $Y=[0,1]^{d}$. The molecular diffusion coefficient $\textbf{A}:Y\to\mathbb{R}^{d\times d}$ is assumed to vary in the cell $Y$, while we consider in this scenario the presence of a volume reaction $\mathcal{R}:\mathbb{R}\to\mathbb{R}$ subject to an internal source $f:\Omega\to\mathbb{R}$. We also take into account the no-flux boundary condition at the internal boundaries, denoted by $\Gamma^{\varepsilon}$, whilst giving the homogeneous Dirichlet boundary condition at the exterior boundary, denoted by $\Gamma^{\text{ext}}$. Essentially, this context can be understood by the following elliptic problem:

[TABLE]

where $L_{\varepsilon}$ is a symmetric operator given by

[TABLE]

1.2 Background and motivation

In this paper, we follow up on our earlier works Khoa and Muntean (2016); Khoa (2017); Khoa et al. (2020) that focus on the asymptotic analysis of semi-linear elliptic problems posed in perforated domains. This well-understood elliptic problem was applied in the studies of the heat transfer in composite materials and of the pressure and phase velocities in porous media flow. Even though those applications are involved in the microscopic system (1), our wish in this paper is to prepare and develop the playground of the model for diffusion, aggregation and surface deposition of a concentration in a porous medium; cf. e.g. Krehel et al. (2015); Timofte (2013). The presence of non-negative scalings stems from our mathematical concerns about the non-trivial asymptotic behaviors of $u_{\varepsilon}$ when $\varepsilon$ tends to 0 and their corresponding rates. In fact, our most recent result in Khoa et al. (2020) has shown that when $\alpha<0$ the macroscopic solution is identically zero after the homogenization process ($\varepsilon\to 0$).

The study of the variable scaled nonlinearities soon appeared in the works Cabarrubias and Donato (2012); Conca et al. (2004), where they considered the homogenization of elliptic problems with a scaled Robin boundary condition. As a concrete motivation of our proposed model (1), the authors in Conca et al. (2004) delved into the chemical reactive flows through the exterior of a domain with distributed reactive obstacles, where the fractional (Langmuir kinetics) and polynomial (Freundlich kinetics) surface reactions were taken into account. Mathematically, such surface reactions posed on $\Gamma^{\varepsilon}$ read as $-\textbf{A}(x/\varepsilon)\nabla u_{\varepsilon}\cdot n=\varepsilon^{\beta}\mathcal{S}(u_{\varepsilon})$ for $\beta\geq 1$, just like what has been studied in Khoa et al. (2020). Note that even though we only consider in (1) the zero Neumann case of the internal boundary $\Gamma^{\varepsilon}$, it is completely similar to adapt our analysis below to the surface reaction (using the same assumption as that of the volume one). Hence, our mathematical analysis in this work is pertinent and applicable.

Arguments obtained from studies of the variable scalings can be helpful in the qualitative analysis of eigen-elements for elliptic boundary value problems with rapidly oscillating coefficients in a perforated cube (cf. e.g. Douanla and Svanstedt (2012); Gryshchuk and de Cristoforis (2013)), while it also concerns low-cost control problems in Muthukumar and Nandakumaran (2009). On top of that, they can be further adapted to complex scenarios that include, for example, the non-stationary Stokes–Nernst–Planck–Poisson system in the dynamics of colloids (cf. e.g. Ray et al. (2012) and references cited therein).

Cf. Khoa et al. (2020), the microscopic solution to the problem $(P_{\varepsilon})$ converges to a macroscopic function that solves a certain non-trivial homogenized problem since the scaling variable $\alpha\geq 0$ is eventually the main factor that determines the presence of the reaction term $\mathcal{R}$ at the macro-scale. In principle, the main purpose of homogenization is to unveil the macroscopic shape of $u_{\varepsilon}$ as it is less time-consuming than any approximation of $u_{\varepsilon}$ itself. Indeed, as is known in the homogenization community, solving $(P_{\varepsilon})$ directly is computationally expensive since the space discretization is inversely proportional to the scale parameter $\varepsilon$. Before the homogenization, one essentially needs the existence and uniqueness of the microscopic solution. Therefore, it ended up with our pursuit to design a linearization scheme that plays these two roles: simultaneously obtaining weak sovability of the microscopic system (for all $\varepsilon>0$) and the linearized homogenized system.

Cf. Khoa and Muntean (2016) as our initiation, a linearization scheme was briefly designed to prove the weak solvability of $(P_{\varepsilon})$ as $\alpha=0$. However, this result was only guaranteed when the diffusion must be very larger than the Lipschitz rate of reactions. Our next evolution in this area went to the work Khoa et al. (2020), where, for the first time, we addressed a linearization scheme for the weak solvability of a semi-linear microscopic system with real variable scalings. Since in those works, our interests were in the asymptotic analysis of the microscopic solution, we did not get through the computational standpoint of the homogenized system. Thereupon, this is the moment to delineate a reliable approximation of $u_{\varepsilon}$ at the macro-scale, which is our second purpose of the aforementioned wish.

1.3 Goals

In this work, we show that the linearization scheme we design is essential in proving the well-posedness of $(P_{\varepsilon})$, but also in deriving the approximate macroscopic solution with certain error estimates. Basically, our theoretical analysis will proceed in accordance with the following diagram:

[TABLE]

Here, $(P_{\varepsilon}^{k})$ denotes the approximate problem of $(P_{\varepsilon})$ by linearization, whilst its macroscopic equation is structured in $(P_{0}^{k})$. The notion behind this approach is to linearize nonlinearities in the model using a suitable choice of the stabilization parameter. As the nonlinearity is supposed to be degenerate at a single point, we are aided by a regularization approach during the linearization process. In this way, we arrive at a regularized form of the nonlinearity, where we can figure out the error estimate between $u_{\varepsilon}$ and $u_{\varepsilon}^{k}$ in $L^{2}$-norm. In our proof, the stabilization is $\varepsilon$-dependent only when the scaling factor $\alpha$ is positive. Meanwhile, $\varepsilon$ does not contribute to the convergence of the linearization scheme for any $\alpha\geq 0$. As $\varepsilon\to 0$ in the homogenization process, the stabilization constant becomes $\varepsilon$-independent for any $\alpha\geq 0$. Henceforth, several stability estimates for the macroscopic solution are easily obtained.

As another vantage of our proposed scheme, we point out that if the solvability of $(P_{\varepsilon})$ and the corresponding macroscopic equation are already known, one can also use the scheme directly to get the approximation $u_{0}^{k}$. In this case, we mean

[TABLE]

Furthermore, in the diagram (10) we can prove that the rate of convergence is $k$-independent (cf. Slodička (2001)), compared to the result we obtain in Corollary 2 for the diagram (6). It is worth mentioning that as we aim to show the unique weak solvability of the microscopic problem by a linearization technique, we follow diagram (6).

1.4 Outline of the paper

The paper is organized as follows. Section 2 is dedicated to introducing notations and necessary assumptions on the input of the problem. In Section 3, we propose an iterations-based variational scheme to linearize the microscopic problem $\left(P_{\varepsilon}\right)$. Accordingly, we obtain the well-posedness of $\left(P_{\varepsilon}\right)$ as well as the rate of convergence by the linearization we choose. Settings of the homogenization are involved in Section 4, where we also state the structures of the cell problems and the limit equation at every step of linearization. Additionally, several types of stability analysis of the scheme at the macro-scale are justified. Section 5 is devoted to the numerical test of the scheme and we close this paper with some concluding remarks in Section 6.

2 Preliminaries

In the sequel, all the constants $C$ are independent of the scaling parameter $\varepsilon$, but their precise values may differ from line to line and may change even within a single chain of estimates. We use either the superscript or subscript $\varepsilon$ to indicate its dependence. Depending on the situation, we denote by $\left|\cdot\right|$ the absolute value of a function, the finite-dimensional Euclidean norm of a vector or the volume of a domain.

The perforated domain $\Omega^{\varepsilon}$ is supposed to be bounded, connected and possesses a $C^{0,1}$-boundary. It is thought to approximate a porous medium and its precise description can be found in Khoa et al. (2015). For brevity, we herein skip the mathematical descriptions of the perforated domain of interest. Instead, we only provide Figure 1 an admissible geometry of our medium and the corresponding microstructure when looking for a graphically schematic representation of the scaling procedure within a natural soil. In addition, the reader can be referred to Hornung and Jäger (1991); Fatima and Muntean (2014); Cancedda (2016) and those collectively mentioned in Section 1 for some concrete results concerning such domains.

Definition 1

(Degenerate class) A real-valued function $F$ is said to be degenerate at a point $x_{0}\in\mathbb{R}$ if we can find $\delta_{1}>0$ such that $0\leq F^{\prime}(x)\leq\delta_{1}$ a.e. in $\mathbb{R}$ and the following conditions hold true:

• •

If $\lim\limits_{x\to s^{+}}F^{\prime}\left(x\right)\lim\limits_{x\to s^{-}}F^{\prime}\left(x\right)=0$ for some $s\in\mathbb{R}$, then $s=x_{0}$.

• •

There exist $\delta_{0},r_{x_{0}}>0$ such that

[TABLE]

where $\mathcal{B}\left(x_{0},\delta_{x_{0}}\right)$ denotes a ball centered at $x_{0}$ with a radius $r_{x_{0}}$.

• •

For $x\in\overline{\mathcal{B}\left(x_{0},r_{x_{0}}\right)}$, $F^{\prime}$ is non-decreasing.

Next, we introduce the space

[TABLE]

equipped with the norm

[TABLE]

Cf. (Cioranescu and Paulin, 1999, Lemma 2.1), one can show the uniform-in-$\varepsilon$ equivalence between this norm and the usual $H^{1}$-norm by the Poincaré-type inequality.

Lemma 1

There exists a constant $C_{p}>0$ independent of $\varepsilon$ such that

[TABLE]

Moreover, we denote by $H^{1}_{0}(\Omega)$ the Hilbert space of weakly differentiable functions $u:\Omega\to\mathbb{R}$ that vanishes on the boundary in the sense of trace. We also use below the space $H^{1}_{\#}(Y_{l})$ to indicate functions in $H^{1}(Y_{l})$ that is $Y$-periodic and has zero mean value.

This is now the moment to state our working assumptions on data involved in the microscopic problem. Those include

$\left(\text{A}_{1}\right)$ The diffusion $\mathbf{A}$ is essentially bounded, $Y$-periodic, symmetric and globally Lipschitz. It satisfies the uniform ellipticity condition in the sense that we can find $\varepsilon$-independent constants $\underline{\gamma},\overline{\gamma}>0$ such that

[TABLE]

$\left(\text{A}_{2}\right)$ The reaction term $\mathcal{R}:\mathbb{R}\to\mathbb{R}$ is degenerate.

$\left(\text{A}_{3}\right)$ The internal source $f$ belongs to $L^{2}(\Omega)$.

Remark 1

The consideration of a degenerate class of the reaction/production term $\mathcal{R}$ starts from our mathematical wish to generalize the previous results in Khoa et al. (2015, 2020), where we exploited linearization schemes under strong assumptions, e.g the global Lipschitz reaction accounting for the Monod equation. In Khoa et al. (2015), we even needed the fact that the fraction of reaction and diffusion bounds has to be sufficiently small. Our generalization in this paper is based upon the molecular interaction of piecewise types. For example, they can be of the power-growth law $F(x)=\left|x\right|^{p},p>1$ in a neighborhood of the degenerate point $x_{0}$ (cf. Definition 1), and follow any global Lipschitz rate outside this neighborhood. A certain example of this typical reaction is introduced in (30).

3 Settings of the iterative variational algorithm

There are several linearization methods investigated in the past and each with its modifications and improvements to serve certain classes of nonlinear partial differential equations. To give a very cursory glance, one may concern the Newton method, cf. Bergamaschi and Putti (1999), whose convergence requires the initial guess to be close to the true solution, albeit its quadratic convergence. The Jäger–Kačur scheme is also renowned for its outstanding performance in the approximations of one-dimensional parabolic problems with a linear convergence; see Jäger and Kačur (1995). In this paper, our method is conventionally in line with the so-called $L$-scheme extensively studied in many distinctive types of parabolic equations; cf. Mitra and Pop (2019) and references cited therein for a short background concerning this typical scheme. In this sense, the so-called stabilization term is added to stabilize the entire linearized equation in the standard variational formulation. Thereby, a linear convergence is obtained under a suitable choice of the stabilization constant.

Definition 2

For each $\varepsilon>0$, a function $u_{\varepsilon}\in H^{1}(\Omega^{\varepsilon})$ is said to be a weak solution to $(P_{\varepsilon})$ if it satisfies

[TABLE]

for all $\varphi\in V^{\varepsilon}$, where $a$: $H^{1}(\Omega^{\varepsilon})\times H^{1}(\Omega^{\varepsilon})\to\mathbb{R}$ is a bilinear mapping given by

[TABLE]

At this stage, we take into account the degeneracy of the reaction term $\mathcal{R}$, especially in the ball where its derivative is zero and non-decreasing (see again in Definition 1). It is worth citing here the Jackson type estimates in the approximation theory of monotone functions by monotone polynomials. In principle, cf. Leviatan (1988), for every monotone non-decreasing function $\bar{f}\in C^{k}[-1,1]$, there are non-decreasing polynomials $p_{n}$, whose degree does not exceed $n$, such that

[TABLE]

where $C$ here is independent of $\bar{f}$ and $n$, and $\omega$ indicates the modulus of continuity of $\bar{f}^{(k)}$. Cf. Kimchi and Leviatan (1976); Roulier (1976) for $\bar{f}\in C^{2k}[-1,1]$ with certain conditions, if $\bar{f}$ possesses strict lower and upper bounds of the available derivatives, then for sufficiently large $n$, the best polynomial approximation to $\bar{f}$ also satisfies the same property. Besides, the corresponding derivatives of the best approximation approach the derivatives of $\bar{f}$, respectively. Cf. Leviatan (2000) for a survey of recent developments of the polynomial approximation, the approximation process can also preserve the monotonicity of $\bar{f}^{\prime}$.

The aforementioned references enable us to assume the existence of a regularization scheme for the degenerate $\mathcal{R}$. In this regard, a function $F_{\gamma}$ for $\gamma>0$ being as a regularization parameter is said to be a regularization of a degenerate function $F$ if one has

• •

$0<\gamma\delta_{0}\leq F_{\gamma}^{\prime}(x)\leq\delta_{1}$ for any $x\in\mathbb{R}$,

• •

$\left|F-F_{\gamma}\right|\leq C\gamma^{\sigma}$ for $\sigma>0$ and for any $x\in\mathbb{R}$.

Technically, the smallness of the lower bound of $F_{\gamma}^{\prime}$ is taken to regularize the degeneracy of $F$ at $x_{0}$ as it is zero at this degenerate point. Thus, one may tacitly look for $F_{\gamma}$ such that $F_{\gamma}^{\prime}(x_{0})=\gamma\delta_{0}$ and should attempt to preserve the “shape” of $F$ through the regularization process. Note that this approach does not mean that the regularization scheme must be a linear mapping; in general, the process can still be nonlinear.

After regularization of the reaction term, the linearized problem in the variational setting is given as follows.

Definition 3

For each $\varepsilon>0$, a linear approximation of $u_{\varepsilon}$ that solves (16) is defined as a sequence $\left\{u_{\varepsilon}^{k}\right\}_{k\in\mathbb{N}^{*}}$ satisfying

[TABLE]

where the stabilization constant $M\geq 0$ and the regularization parameter $\gamma_{k}>0$ are selected later. The initial guess is taken as $u_{\varepsilon}^{0}=0$.

In the following theorem, we prove that this sequence is well-defined and exists uniquely in $H^{1}(\Omega^{\varepsilon})$.

Theorem 3.1

Assume $\left(\text{A}_{1}\right)$–$\left(\text{A}_{3}\right)$ hold. By choosing $M=\eta+\varepsilon^{\alpha}\delta_{1}$ for $\eta>0$, the approximate problem (19) admits at most a weak solution in $H^{1}(\Omega^{\varepsilon})$.

Proof

It suffices to consider the first-loop problem of (19) (i.e. $k=1$), which reads as

[TABLE]

where $\tilde{f}:=f-\varepsilon^{\alpha}\mathcal{R}_{\gamma_{k}}(0)$. Hereby, we can introduce the bilinear form $\mathcal{B}:H^{1}\left(\Omega^{\varepsilon}\right)\times H^{1}\left(\Omega^{\varepsilon}\right)\to\mathbb{R}$ given by $\mathcal{B}(u,\varphi):=a(u,\varphi)+M\left\langle u,\varphi\right\rangle_{L^{2}\left(\Omega^{\varepsilon}\right)}$. Hence, we complete the proof of the theorem by virtue of the natural $\varepsilon$-independent coerciveness and continuity of $\mathcal{B}$ in $H^{1}(\Omega^{\varepsilon})$.

Lemma 2

Let $\left\{p_{k}\right\}_{k\in\mathbb{N}^{*}}$ and $\left\{q_{k}\right\}_{k\in\mathbb{N}^{*}}$ be sequences of nonnegative real numbers that obey the following recursion

[TABLE]

where $a_{k}$ and $b_{k}$ are also nonnegative real numbers. Then, it holds

[TABLE]

Proof

The proof is trivial and it can be found in A.

Theorem 3.2

Let $\alpha\geq 0$. Assume $\left(\text{A}_{1}\right)$–$\left(\text{A}_{3}\right)$ hold and let $\gamma_{k}>0$ be non-increasing with $\sigma>0$. Then by choosing $M=\eta+\varepsilon^{\alpha}\delta_{1}$ for $\eta>0$, the sequence $\left\{u_{\varepsilon}^{k}\right\}_{k\in\mathbb{N}^{*}}$ of the variational problem (19) possesses the following property

[TABLE]

Proof

Consider $w_{\varepsilon}^{k}:=u_{\varepsilon}^{k}-u_{\varepsilon}^{k-1}$ for $k\geq 2$ as a difference function between the $k$th and $(k-1)$th steps of approximation. Then the difference equation is provided by

[TABLE]

Define $g_{\gamma_{k}}\left(t\right):=\varepsilon^{\alpha}\mathcal{R}_{\gamma_{k}}\left(t\right)-Mt$ for $t\in\mathbb{R}$. By taking the test function $\varphi=w_{\varepsilon}^{k}$, (24) becomes

[TABLE]

Observe that $\left|g^{\prime}_{\gamma_{k}}\right|\leq M-\varepsilon^{\alpha}\gamma_{k}\delta_{0}$ and in view of the fact that

[TABLE]

resulting from the regularization factor we apply, we estimate that

[TABLE]

Upon the monotonicity of $\gamma_{k}$, we use the Young inequality to get

[TABLE]

Using again the Young inequality, we also obtain the upper bound of $I_{2}$ as follows:

[TABLE]

Thereby, after some rearrangements, we find that

[TABLE]

At present, we take $M=\eta+\varepsilon^{\alpha}\delta_{1}$ for $\eta>0$ (independent of $\varepsilon$ and $k$) in (26). Then we apply the Poincaré inequality (cf. Lemma 1) to arrive at

[TABLE]

and it thus follows from Lemma 2 that

[TABLE]

where we have denoted by

[TABLE]

Naturally, $b_{k}\in(0,1)$ and thus the product of $b_{i}$ approaches 0 when $k$ tends to infinity in the sense that

[TABLE]

This completes the proof of the theorem.

In (23), we observe that the stability of the scheme is essentially dependent on the partial sums of the series of $a_{k}$, particularly of the choice of the regularization parameter $\gamma_{k}$. Since in this paper we obtain a strong convergence along with a particular error estimate, harmonic series are not reliable in ensuring that the sequence $\left\{u_{\varepsilon}^{k}\right\}_{k\in\mathbb{N}^{*}}$ is Cauchy in $H^{1}\left(\Omega^{\varepsilon}\right)$. This is hindered by the convergence-towards-zero of the series of $a_{k}$ when using the standard triangle inequality. In simpler terms, harmonic series (or even hyperharmonic series) are mostly either divergent or convergent to a non-zero constant. It is worth noting that the product of $b_{k}$ possesses an exponential-like decay as $k\to\infty$. Then the same behavior should be applied to the series of $a_{k}$ by looking for a geometric progression of $a_{k}$.

As a concrete example, we take into account the power-law reaction rate in a unit domain, which reads as

[TABLE]

Therefore, in the interval $[0,1]$ we can choose a regularization of $\mathcal{R}$ as follows:

[TABLE]

leading to the fact that $\frac{\delta_{0}}{2^{(p-1)^{2}+(k+1)(p-1)}}\leq\mathcal{R}_{\gamma_{k}}^{\prime}\leq 1$ and the following estimate

[TABLE]

In this way, we indicate that $\gamma_{k}=\frac{1}{2^{(p-1)^{2}+(k+1)(p-1)}}$ and $\sigma=1+\frac{1}{p-1}>1$. It is worth noting that

[TABLE]

and thereupon, we, in accordance with the estimate (23), provide the following theorem.

Theorem 3.3

Under the assumptions of Theorem 3.2, if we further choose $\gamma_{k}$ in such a way that

[TABLE]

and suppose that

[TABLE]

Then, the iterative sequence $\left\{u_{\varepsilon}^{k}\right\}_{k\in\mathbb{N}^{*}}$ is Cauchy in $H^{1}\left(\Omega^{\varepsilon}\right)$. Consequently, the microscopic problem $(P_{\varepsilon})$ admits a unique solution $u_{\varepsilon}$ in $H^{1}\left(\Omega^{\varepsilon}\right)$.

Proof

According to Theorem 3.2, it is straightforward to find an $\varepsilon$-independent upper bound for (23) as follows:

[TABLE]

where we have essentially used the binomial identity. Thereby, for any $k,r\in\mathbb{N}^{*}$ we obtain

[TABLE]

using the standard triangle inequality. This way we prove that $\left\{u_{\varepsilon}^{k}\right\}_{k\in\mathbb{N}^{*}}$ is Cauchy in $H^{1}\left(\Omega^{\varepsilon}\right)$. Consequently, there exists a unique function $u_{\varepsilon}\in H^{1}\left(\Omega^{\varepsilon}\right)$ such that $u_{\varepsilon}^{k}\to u_{\varepsilon}$ as $k\to\infty$. Furthermore, it is straightforward to get that

[TABLE]

Thus, it enables us to confirm the existence and uniqueness of $u_{\varepsilon}$ to the microscopic problem $(P_{\varepsilon})$ in $H^{1}\left(\Omega^{\varepsilon}\right)$. This completes the proof of the theorem.

Corollary 1

Under the assumptions of Theorem 3.3, the following rate of convergence holds

[TABLE]

Proof

The proof of this corollary is obvious by the aid of (36) when taking $r\to\infty$.

Remark 2

It is worth mentioning here that if we already know that $u_{\varepsilon}\in H^{1}(\Omega^{\varepsilon})$ is a unique solution of the microscopic problem $(P_{\varepsilon})$, then the partial sums of the series $\sum_{j=2}^{k-1}a_{j}\prod_{i=j+1}^{k}b_{i}$ (cf. (28)) just needs to approach 0 as $k\to\infty$ with a certain convergence rate, by considering the energy-like estimate of the difference between the linearized problem and $(P_{\varepsilon})$. Accordingly, the whole error estimate is controlled by such a rate of that partial sums. This exactly mimics the proof in Slodička (2001) where a minimal polynomial rate $k^{-\omega}$ for $\omega\in(0,1)$, which basically leads to the harmonic progression, is sufficiently taken into account.

Additionally, it is straightforward to obtain the stability analysis of the scheme $\left\{u_{\varepsilon}^{k}\right\}_{k\in\mathbb{N}^{*}}$ we construct above, which, in principle, provides concretely $\varepsilon$-independent a priori estimates in $H^{1}(\Omega^{\varepsilon})$. This somewhat enables us to get the existence of a weak solution to the problem $(P_{\varepsilon})$ in $H^{1}(\Omega^{\varepsilon})$ by the standard compactness argument, if we are able to derive, at least, the weak convergence of the reaction term in $L^{2}(\Omega^{\varepsilon})$ after passing to the limit. However, the uniqueness result may not always be achievable by this strategy. In this work, we do not go beyond this matter and will leave it for the future works.

4 Settings of the homogenization

In the previous section, the rigorous error estimate for the linearization scheme has been obtained in $H^{1}(\Omega^{\varepsilon})$; cf. Corollary 1. In this section, we only exploit the $L^{2}$ error estimate, although it is well-known from the corrector estimate for the homogenization limit that the microscopic solution of a linear elliptic equation approaches the macroscopic solution of the corresponding homogenized elliptic equation with a rate $\mathcal{O}\left(\varepsilon^{\frac{1}{2}}\right)$ in $H^{1}(\Omega^{\varepsilon})$; cf. (Cioranescu and Paulin, 1999, Corollary 2.29). From here on, it is very easy to adapt this result since our approximate problem defined in (3) is all linear at every step $k$.

Our goal here is to introduce the structures of the non-trivial homogenized problem $(P_{0})$ for $(P_{\varepsilon})$ as well as its cell problems for the sake of computations in Section 5. When doing so, we remark that $u_{\varepsilon}$ satisfies an a priori estimate by means of $\left\|u_{\varepsilon}\right\|_{H^{1}\left(\Omega^{\varepsilon}\right)}\leq C$ established in Theorem 3.1 and by taking into account the usual zero extension on $u_{\varepsilon}$ from $H^{1}(\Omega^{\varepsilon})$ to $H^{1}(\Omega)$. Accordingly, we only need to take the test function $\varphi=\psi_{0}(x)+\varepsilon\psi_{1}\left(x,\frac{x}{\varepsilon}\right)$ for $\psi_{0},\psi_{1}\in\mathcal{D}(\Omega;C^{\infty}_{\#}(Y))$ in (3). Henceforward, the compactness result allows us to extract subsequences from bounded sequences and to obtain the passage to the two-scale limit. For detailed results concerning the two-scale convergence method for the linear elliptic equation, we refer the interested reader to e.g. Henning and Ohlberger (2009); Ray et al. (2012) under the theoretical results of the two-scale convergence postulated in Allaire (1992); Nguetseng (1989).

After plugging that typical test function we have

[TABLE]

then by passing to the limit $\varepsilon\to 0$, we are led to the following limit equation:

• •

Case $\alpha>0$:

[TABLE]

• •

Case $\alpha=0$:

[TABLE]

These two cases are almost the same since we are at the linearization stage. Typically, we choose $\psi_{0}=0$ in both two cases to get

[TABLE]

where we have used the integration by parts with respect to $y$. In this way, we use the separation of variables to find

[TABLE]

Here, $\chi^{k}_{i}$ is called as the cell function that solves the following cell problem:

[TABLE]

It is worth noting that the cell problem for every step $k$ remains unchanged, so that our computations will be less expensive in the sense that we do not need to compute the vector-valued $\chi=\chi^{k}=\left(\chi_{i}^{k}\right)_{1\leq i\leq d}$ at every $k$. Furthermore, one can prove that such $\chi^{k}\in H^{1}_{\#}(Y_{l})$. Due to the non-convexity of $Y_{l}$, the regularity of the unique function $\chi$ stops at $H^{1+s}(Y_{l})$ for $s\in(-1/2,1/2)$, no matter how smooth the involved data are; see Savaré (1998) for detailed concerns.

Now, choosing $\psi_{1}=0$ and then applying the integration by parts with respect to $x$, we obtain the equation for $u_{0}^{k}$ as follows:

[TABLE]

posed in $\Omega$. Here, this equation is endowed with the Dirichlet boundary condition $u_{0}^{k}=0$ at $\partial\Omega$ and $A^{0}$ is known as the homogenized coefficient given by

[TABLE]

where $\delta_{ij}$ stands for the constants of the identity matrix. Additionally, this coefficient satisfies the uniform ellipticity condition by virtue of the well-known Voigt–Reiss inequality; see e.g. Sviercoski et al. (2008). Thus, proof of the well-posedness of the macroscopic problem in $H^{1}_{0}(\Omega)$ is standard.

Remark 3

It is worth noting that the linearization scheme of the microscopic system (1) plays two roles in this paper. As the primary goal, it allows us to prove the well-posedness of the microscopic system for any $\varepsilon>0$. The second purpose of this scheme is to obtain the corresponding approximate homogenized system (41). To this end, this finding enables us to compute the reliable homogenized solution. This is a new story of the linearization scheme developed in Slodička (2001) for semi-linear problems at only the macro-scale, where the author proved its strong convergence by a priori knowledge of the well-posedness of the original problem. Moreover, this is a novel contribution of the recent result Khoa et al. (2020), where the authors only studied the global Lipschitz reaction and the asymptotic analysis of (1).

Remark 4

As we now work on the linearized elliptic problem at the macro-scale, we below follow the same strategy of Theorems 3.2, 3.3 and those of Slodička (2001) to guarantee the stability (in some suitable norms) of the scheme in the macroscopic scaling and with respect to the number of iteration $k$.

Theorem 4.1

Under the assumptions of Theorem 3.3 and suppose that

[TABLE]

where $c_{p}>0$ is the standard Poincaré111For any $u\in H_{0}^{1}\left(\Omega\right)$, it holds $\left\|u\right\|^{2}_{L^{2}\left(\Omega\right)}\leq c_{p}\left\|\nabla u\right\|^{2}_{L^{2}\left(\Omega\right)}$. constant. Then, the iterative sequence $\left\{u_{0}^{k}\right\}_{k\in\mathbb{N}^{*}}$ is Cauchy in $H^{1}_{0}(\Omega)$.

Proof

It is trivial to prove that the functions $u_{0}^{k}$ are all in $H_{0}^{1}(\Omega)$. Observe that (41) is structured by the cases of $\alpha$, we thus also divide the proof here into two parts. In the first part, we treat the following equation:

[TABLE]

Recall that the problem under consideration is associated with the zero Dirichlet boundary condition.

We now use the test function $\bar{\varphi}\in H^{1}_{0}(\Omega)$ to arrive at the following variational formulation:

[TABLE]

where $\bar{a}:H_{0}^{1}(\Omega)\times H_{0}^{1}(\Omega)\to\mathbb{R}$ is a bilinear mapping given by

[TABLE]

Then, it is straightforward to compute the difference equation by putting $v_{0}^{k}=u_{0}^{k}-u_{0}^{k-1}$. This function essentially satisfies the following equation:

[TABLE]

and then, by taking $\bar{\varphi}=v_{0}^{k}$ it leads to

[TABLE]

Due to the standard Poincaré inequality, we have

[TABLE]

and by using Lemma 2, we estimate that

[TABLE]

Similar to proof of Theorem 3.3, we can prove that $\left\{u_{0}^{k}\right\}_{k\in\mathbb{N}^{*}}$ is Cauchy in $H^{1}_{0}(\Omega)$.

Next, we consider the following equation:

[TABLE]

We proceed as above by taking into account the variational formulation of the difference equation. By so doing, (4) becomes

[TABLE]

where we have denoted by $\bar{g}_{\gamma_{k}}(t)=\mathcal{R}_{\gamma_{k}}(t)-\left(\eta+\delta_{1}\right)t$.

Since $\left|\bar{g}_{\gamma_{k}}^{\prime}\right|\leq\eta+\delta_{1}-\gamma_{k}\delta_{0}$ and $\left|\mathcal{R}_{\gamma_{k-1}}-\mathcal{R}_{\gamma_{k}}\right|\leq C\gamma_{k-1}^{\sigma}$, we then apply the Young inequality to obtain

[TABLE]

and after some rearrangements and applying the standard Poincaré inequality, we arrive at

[TABLE]

Henceforward, with the aid of Lemma 2 we have the same estimate as (23) and by the choice of $\gamma_{k}$ in Theorem 3.3, we consequently prove that $\left\{u_{0}^{k}\right\}_{k\in\mathbb{N}^{*}}$ is Cauchy in $H^{1}_{0}(\Omega)$. This completes the proof of the theorem.

As argued in Remark 2, the geometric progression is required to prove the existence and uniqueness of the microscopic problem $\left(P_{\varepsilon}\right)$ in $H^{1}(\Omega^{\varepsilon})$, and to avoid the case $\alpha>0$ where the error bound is arbitrarily slow. Since we are now dealing with the macroscopic framework, we can get a better “stability analysis”. This argument is shown in the following theorem as a stability analysis of the macroscopic scheme in $L^{\infty}(\Omega)$. Additionally, this can be applied to Theorem 4.1 where the assumption (43) is no longer necessary by not using the Poincaré inequality.

Theorem 4.2

Consider the case $\alpha=0$. Suppose that the internal source $f$ and the regularization $\mathcal{R}_{\gamma_{k}}$ are smooth. Under the assumptions of Theorem 4.1, if we can choose

[TABLE]

the iterative sequence $\left\{u_{0}^{k}\right\}_{k\in\mathbb{N}^{*}}$ is stable in the sense that

[TABLE]

Proof

In this proof, we also follow the same vein as proof of Theorem 4.1 where we consider two typical structures of the limit equation (41) and investigate their corresponding difference equation by the standard variational setting. For brevity, we do not mention again those equations, but the governing difference equations. Due to the linear problem as well as the smoothness of $f$ and $\mathcal{R}_{\gamma_{k}}$, it is also trivial to prove that the functions $u_{0}^{k}$ are in $L^{\infty}(\Omega)$, cf. Agmon et al. (1959). In this regard, we recall the difference equation

[TABLE]

where the test function $\bar{\varphi}\in H_{0}^{1}(\Omega)\cap L^{\infty}(\Omega)$ is now taken into account.

Now we set

[TABLE]

and put

[TABLE]

Then we assume that $\left|W_{k}^{1}\right|>0$. By choosing $\bar{\varphi}=\left(v_{0}^{k}+D_{k}\right)^{-}$ in (53) where $f^{-}:=\min\left\{f,0\right\}$, we project (53) from $\Omega$ to the set $W_{k}^{1}$. Thus, below we will accompany the notation $W_{k}^{1}$ to indicate the fact that we are working in that set, although the essential integral has to be posed in $\Omega$. In fact, (53) now becomes

[TABLE]

At this stage, we see that the bilinear form $\bar{a}_{W_{k}^{1}}\left(v_{0}^{k},\bar{\varphi}\right)$ is non-negative certainly due to the presence of the gradient. Furthermore, it holds for a.e. $x\in W_{k}^{1}$ that

[TABLE]

which implies

[TABLE]

This essentially contradicts the equality (56) and thus the assumption $\left|W_{k}^{1}\right|>0$ does not hold. Equivalently, it means that for a.e. $x\in\Omega$, it holds $u_{0}^{k}-u_{0}^{k-1}\geq-D_{k}$.

In the same vein, we can prove that $u_{0}^{k}-u_{0}^{k-1}\leq D_{k}$ by setting

[TABLE]

and employing the test function $\bar{\varphi}=\left(v_{0}^{k}-D_{k}\right)^{+}$ where $f^{+}:=\max\left\{f,0\right\}$. For brevity, we omit the details and leave it to the reader.

Henceforward, we obtain

[TABLE]

In view of the fact that $\left|\bar{g}_{\gamma_{k}}^{\prime}\right|\leq\eta+\delta_{1}-\gamma_{k}\delta_{0}$ and $\left|\mathcal{R}_{\gamma_{k-1}}-\mathcal{R}_{\gamma_{k}}\right|\leq C\gamma_{k-1}^{\sigma}$, it reveals

[TABLE]

Thanks to Lemma 2, one then deduces

[TABLE]

where

[TABLE]

It suffices to take into account the second term on the right-hand side of (61). Observe that we can bound it from above by using the standard inequality $1+x\leq\text{exp}(x)$ for $x\in\mathbb{R}$. Indeed, it is aided by the sum and integral inequality222$\int_{j}^{k}\gamma\left(x\right)dx\leq\sum_{i=j+1}^{k}\gamma\left(i\right)\leq\int_{j+1}^{k+1}\gamma\left(x\right)dx$ for any non-decreasing $\gamma$. that

[TABLE]

We now employ the elementary inequality $\left(j+2\right)^{\sigma}\leq 2^{\sigma-1}\left(j^{\sigma}+2^{\sigma}\right)\leq Cj^{\sigma}$ for $j\geq 2$ to arrive at

[TABLE]

by virtue of $\sigma>\frac{\delta_{0}}{\eta+\delta_{1}}$ and $1<\frac{k+2}{k+1}<2$.

Hence, this way we complete the proof of the theorem.

Since the homogenized system is linearized and associated with constant coefficients, one can obtain the error estimate between $u_{\varepsilon}^{k}$ and $u_{0}^{k}$ in $H^{1}(\Omega^{\varepsilon})$ by adapting the classical two-scale asymptotic expansion

[TABLE]

Theorem 4.3

Assume $\left(\text{A}_{1}\right)$ holds. Suppose that $\partial\Omega\in C^{6}$ and $f,\mathcal{R}(0)\in H^{4}(\Omega)$. Then there exists $C>0$ independent of $\varepsilon$ such that

[TABLE]

Proof

The proof was already obtained via (Cioranescu and Paulin, 1999, Corollary 2.30). Here, we single out how that proof works in our particular case. Starting from $k=1$, our solution $u_{0}^{1}$ belongs to $H^{6}(\Omega)$ using the standard higher boundary regularity result for elliptic equations (cf. (Evans, 2010, Section 6.3)). Since $\mathcal{R}\left(0\right)\in H^{4}\left(\Omega\right)$, we get $\mathcal{R}_{\gamma_{k}}\left(u_{0}^{k-1}\right)\in H^{4}(\Omega)$ in the sense that

[TABLE]

Therefore, we conclude that any $u_{0}^{k}$ belongs to $H^{6}(\Omega)$. The embedding $H^{6}(\Omega)\subset C^{4}(\overline{\Omega})$ guarantees the fact that the derivatives of $u_{0}^{k}$ up to the fourth order are bounded. Moreover, recall that cf. Savaré (1998), the cell functions $\chi^{k}$ are in $H^{1}(Y_{l})$. Hence, we can adapt the proof of (Cioranescu and Paulin, 1999, Corollary 2.30) to prove that (63) holds true.

Ultimately, we deduce the following error estimate between $u_{\varepsilon}$ and $u_{0}^{k}$.

Corollary 2

Under the assumptions of Theorem 3.3 and Theorem 4.3, the following rate of convergence holds

[TABLE]

where $C_{k}=C\left(\left\|u_{0}^{k-1}\right\|_{H^{4}\left(\Omega\right)}\right)>0$.

In Corollary 2, we see that due to the presence of $C_{k}$ it is not necessary to take $k$ very large to get a fine approximation of $u_{\varepsilon}$. Theoretically, this $C_{k}$ largeness might slow down the smallness of the quantity $\varepsilon^{1/2}$. From the computational standpoint, usually one should stop at $k=5$ at most because the algorithm is robust and it avoids being time-consuming.

5 Numerical implementation

In this section, we numerically investigate the potential of the developed iterative method in approximating nonlinear elliptic problems, which are described on complex porous domains. Here, we focus on the case $\alpha=0$ since this case not only remains nonlinear in the macroscopic problem, but also can be considered as a paradigm for treating other problems of interest. Besides, the case, $\alpha>0$, is simple and can be handled in a more straightforward manner since the reaction term converges to [math] as $\varepsilon$ tends [math].

For the purpose of illustration, let us consider the problem (1) in a two-dimensional unit square $\Omega=(0,1)\times(0,1)$. The highly oscillatory diffusion coefficient is chosen as

[TABLE]

and we take $\mathcal{R}(u)$ as in (30) with $p=2$, and the source term is $f=1$.

In order to demonstrate the $L^{2}$ estimate of the error between $u_{\varepsilon}$ and $u_{\varepsilon}^{k}$, we first solve the problem (41) for $u_{\varepsilon}^{k}$, which involve solving (40) and (42) for the cell functions $\chi_{i},i=1,2,$ and for the effective diffusion coefficient, $A^{0}$, respectively. We consider (40) in the unit cell $Y=(0,1)\times(0,1)$ with a hole of radius $r=0.4$ and porosity, $|Y_{l}|=1-\pi r^{2}$. Moreover, we take the constants $\eta=0.4$ and $\delta_{1}=1$ for the stabilization constant $M$, cf. Theorem 3.1. Note that the regularization $\mathcal{R}$ is given, according to (30), by

[TABLE]

Eventually, by plugging the cell solutions $\chi_{i},i=1,2$ in (42), we can compute the homogenized diffusion coefficient, as follows:

[TABLE]

The original problem (1) for $u_{\varepsilon}$ is solved using the Newton–Raphson method and the P1 standard finite elements on a non-uniform mesh discretization with a size $h$, satisfying $h<\varepsilon$. Since the mesh is non-uniform, we denote by $\max_{h}$ the largest value of the mesh size to qualify the condition $h<\varepsilon$. As $\varepsilon$ decreases from 0.5, $\max_{h}$ also decreases and remains consistent with the requirement, $h<\varepsilon$, for multiscale simulation of homogenization problems as shown in Table 1. The iterative scheme for problems (19) and (41) are solved until the difference between the $L^{2}$ estimates of successive iterations is close to zero.

To assess the efficiency of the linearization of (1) in (19), we take into account the following error estimates:

[TABLE]

The relative error between $u_{\varepsilon,h}^{k}$ and $u_{0,k}^{h}$ are tabulated in Table 2. It indicates that the discrepancy between the solutions of the linearized microscopic problem and the linearized macroscopic problem is good within $k=4$. However, the convergence rate of the Newton’s iteration for $u_{\varepsilon,h}$ is fast when compared with the convergence rate of the linearization algorithm for $u^{k}_{\varepsilon,h}$. One key result, as demonstrated in Tables 1 and 2, is that the proposed linearization algorithm conveniently promotes the passage from the nonlinear microscopic description to the corresponding macroscopic description, without the need of performing Taylor’s expansion for nonlinear terms as is customary in the classical homogenization theory. We also remark that the Newton’s iteration usually needs a fine initial guess to attain convergence, while the choice is arbitrary for our linearization scheme.

In addition, in Table 1, the norm of the difference between $\nabla u_{\varepsilon,h}$ and $\nabla u^{k}_{0,h}$ in $L^{2}$ are listed. When $\varepsilon$ is decreasing, the error becomes smaller with $k$ up to 2. It confirms that the performance of our technique on the problem is good enough for such kind of classical homogenization problems. However, a stronger $H^{1}$ error estimate could be observed by considering terms of correction in the estimate, which is beyond the scope of the present manuscript.

As to our numerical results, by fixing $\varepsilon=0.1$ the relative $L^{2}$ error between $u_{\varepsilon,h}$ and $u_{\varepsilon,h}^{k}$ with $k=2$ gives $1.743\times 10^{-2}$ and we have illustrated these solutions in Figure 2. Figure 3 depicts our linearized macroscopic solution at $k=4$ and the microscopic solution when $\varepsilon$ varies from $0.25$ to $0.025$. We can conclude from Figures 2 and 3 that our linearization scheme performs very well in both microscopic and macroscopic contexts. Moreover, the decrease in the error when $\varepsilon$ becomes smaller, as demonstrated in Tables 1 and 2, shows the consistency of our method with the classical homogenization theory, i.e. the convergence of $u_{\varepsilon}$ to $u_{0}$, for some $k$ in the iterative scheme.

6 Concluding remarks

In this work, we have proposed a regularization- and linearization-based scheme to construct efficient approximations of both microscopic and macroscopic problems. Although, for example, in power-law nonlinearity, a geometric regularization parameter is needed to prove the well-posedness of microscopic problem, in practice one can utilize the harmonic progression at the macro-scale to get convergence from the scheme without paying attention to the rate. Note that the harmonic choice of the regularization parameter is rather well-suited to polynomial approximations as deduced in (18). We emphasize that the arguments in this paper can be typically applied to single out a reliable approximation of the macroscopic equation if the weak solvability of the microscopic scenario is mathematically known. This approach could be helpful for engineering needs, among several types of linearization methods. Furthermore, in upcoming works we will attempt to adopt the so-called boundary layers correctors (see Versieux and Sarkis (2006)) to our context to improve the error estimates in the numerical perspective.

Appendix A Auxiliary proofs

Proof of Lemma 2

The proof is essentially done by induction. Indeed, it trivially holds for $k=3$. For any $k=n$, one now suppose that

[TABLE]

Our aim is to prove that it still holds true for the case $k=n+1$, i.e.

[TABLE]

Using (67), we derive that

[TABLE]

which completes the proof of the lemma.

Acknowledgements

V.A.K. thanks Prof. Adrian Muntean for being his supervisor since March, 2015 and for giving him invaluable advice. V.A.K thanks Prof. Iuliu Sorin Pop (Hasselt, Belgium) for recent supports in his research career and acknowledges the hospitality of the Hasselt University during the time he is hosted as a postdoctoral fellow. N.N.N. acknowledges the support of the project INdAM Doctoral Programme in Mathematics and/or Applications Cofunded by Marie Sklodowska-Curie Actions, acronym: INdAM-DP-COFUND-2015, grant number: 713485.