Rate of Convergence of General Phase Field Equations towards their Homogenized Limit

TL;DR

This paper derives the first error estimates for fourth order, nonlinear phase-field equations in heterogeneous porous media, revealing a convergence rate of approximately     in terms of heterogeneity scale.

Contribution

It provides the first derivation of error estimates for fourth order, homogenized, nonlinear phase-field equations in strongly heterogeneous materials.

Findings

Error estimate with convergence rate     for the homogenized solution

First analysis of fourth order phase-field equations in heterogeneous media

Motivates new approaches in modeling and computation of interfacial phenomena

Abstract

Over the last few decades, phase-field equations have found increasing applicability in a wide range of mathematical-scientific fields (e.g. geometric PDEs and mean curvature flow, materials science for the study of phase transitions) but also engineering ones (e.g. as a computational tool in chemical engineering for interfacial flow studies). Here, we focus on phase-field equations in strongly heterogeneous materials with perforations such as porous media. To the best of our knowledge, we provide the first derivation of error estimates for fourth order, homogenized, and nonlinear evolution equations. Our fourth order problem induces a slightly lower convergence rate, i.e., $\epsilon^{1/4}$, where $\epsilon$ denotes the material's specific heterogeneity, than established for second-order elliptic problems (e.g. \cite{Zhikov2006}) for the error between the effective macroscopic solution of the (new) upscaled formulation and the solution of the microscopic phase field problem. We hope that our study will motivate new modelling, analytic, and computational perspectives for interfacial transport and phase transformations in strongly heterogeneous environments.

Equations

$$\begin{array}{r}e(\varphi ):=\frac{1}{\lambda }F(\varphi )+\frac{\lambda }{2}{\left|\nabla \varphi \right|}^2,\end{array}$$

$$\begin{split}{\color[rgb]{0,0,0}F(u)=\sum_{i=2}^{2r}b_{i}u^{i}\,,\quad ib_{i}=a_{i-1}\,,\quad 2\leq i\leq 2r}\,,\end{split}$$

$$\begin{split}{\color[rgb]{0,0,0}f(u)=\sum_{i=1}^{2r-1}a_{i}u^{i}\,,\quad r\in\mathbb{N}\,,\quad r\geq 2}\,,\end{split}$$

$$g(\mathbf{x})=-\frac{\gamma }{C_h}a(\mathbf{x}),$$

$$\begin{array}{r}\text{(Homogeneous case)}\{\begin{array}{ll}\frac{\partial }{\partial t}\varphi =\mathrm{div}\left(\hat{\mathrm{M}}\nabla \left(\frac{1}{\lambda }f(\varphi )-\lambda \Delta \varphi \right)\right)& \text{in }{\Omega }_T,\\ {\nabla }_n\varphi :=\mathbf{n}\cdot \nabla \varphi =g(\mathbf{x})& \text{on }\partial {\Omega }_T,\\ {\nabla }_n\Delta \varphi =0& \text{on }\partial {\Omega }_T,\end{array}\end{array}$$

$$\begin{array}{r}\mu (\varphi ):=-\lambda \Delta \varphi +\frac{1}{\lambda }f(\varphi ),\end{array}$$

$$\begin{array}{r}\text{Hele-Shaw/Mullins-Sekerka problem:}\{\begin{array}{ll}\Delta \mu =0& \text{in }\Omega \setminus {\Gamma }_t,\\ \mathbf{n}\cdot \nabla \mu =0& \text{on }\partial \Omega ,\\ \mu =\sigma \kappa & \text{on }{\Gamma }_t,\\ v=\frac{1}{2}{\left[\mathbf{n}\cdot \nabla \mu \right]}_{{\Gamma }_t}& \text{on }{\Gamma }_t,\\ {\Gamma }_0={\Gamma }_{00}& \text{if }t=0,\end{array}\end{array}$$

$$\begin{array}{rl}{\Omega }^{ϵ}& :=\bigcup _{\mathbf{z}\in {\mathbb{Z}}^d}ϵ\left(Y^1+\mathbf{z}\right)\cap \Omega ,B^{ϵ}:=\bigcup _{\mathbf{z}\in {\mathbb{Z}}^d}ϵ\left(Y^2+\mathbf{z}\right)\cap \Omega =\Omega \setminus {\Omega }^{ϵ},\end{array}$$

$$\begin{array}{r}\text{(Micro porous case)}\{\begin{array}{ll}{\partial }_t{\varphi }_{ϵ}=\mathrm{div}\left(\hat{\mathrm{M}}\nabla \left(-\lambda \Delta {\varphi }_{ϵ}+\frac{1}{\lambda }f({\varphi }_{ϵ})\right)\right)& \text{in }{\Omega }_T^{ϵ},\\ {\nabla }_n{\varphi }_{ϵ}:=\mathbf{n}\cdot \nabla {\varphi }_{ϵ}=0& \text{on }\partial {\Omega }_T^{ϵ},\\ {\nabla }_n\Delta {\varphi }_{ϵ}=0& \text{on }\partial {\Omega }_T^{ϵ},\\ {\varphi }_{ϵ}(\mathbf{x},0)=\psi (\mathbf{x})& \text{on }{\Omega }^{ϵ}.\end{array}\end{array}$$

$$H^{2}_{E}(\Omega):=\left\{\phi\in H^{k}(\Omega)\,\bigl{|}\,\nabla_{n}\phi=0\text{ and }\overline{\phi}{\color[rgb]{0,0,0}:=\frac{1}{|\Omega|}\int_{\Omega}\phi\,d{\bf x}}=0,\,{\color[rgb]{0,0,0}k\geq 2}\right\}$$

$$\left(-\Delta \varphi ,\phi \right)=\left(-\Delta (-\Delta {)}^{-1}w,\phi \right)=\left(w,\phi \right),$$

$$\begin{array}{rl}\left({\partial }_t(-\Delta {)}^{-1}w,\phi \right)-\left(\lambda \mathrm{div}\left(\hat{\mathrm{M}}\nabla w\right),\phi \right)& =\left(\mathrm{div}\left(\frac{\hat{\mathrm{M}}}{\lambda }\nabla f(\varphi )\right),\phi \right),\\ \left(\nabla \varphi ,\nabla \phi \right)& =\left(w,\phi \right),\end{array}$$

$$\begin{array}{r}\text{(Splitting)}\{\begin{array}{ll}{\partial }_t(-\Delta {)}^{-1}w-\lambda \mathrm{div}\left(\hat{\mathrm{M}}\nabla w\right)=\mathrm{div}\left(\frac{\hat{\mathrm{M}}}{\lambda }\nabla f(\varphi )\right)& \text{in }{\Omega }_T,\\ {\nabla }_{n}w=-{\nabla }_n\Delta \varphi =0& \text{on }\partial {\Omega }_T,\\ -\Delta \varphi =w& \text{in }{\Omega }_T,\\ {\nabla }_n\varphi =0& \text{on }\partial {\Omega }_T,\\ \varphi (\mathbf{x},0)=\psi (\mathbf{x})& \text{in }\Omega .\end{array}\end{array}$$

$$\begin{array}{r}{k}_0{\left|u\right|}^{\alpha -2}-k_1\le f^{\prime }(u)\le k_2{\left|u\right|}^{\alpha -2}+k_3.\end{array}$$

$$\begin{split}\left({\color[rgb]{0,0,0}f(a)-f(b)},a-b\right)&\geq a_{1}\left({\color[rgb]{0,0,0}f^{\prime}(a)}(a-b),a-b\right)-a_{2}\left|a-b\right|^{2+a_{3}}\quad\forall\left|a\right|\leq 2a_{4}\,,\\ aF^{\prime\prime}(a)&\geq 0\qquad\forall\left|a\right|\geq a_{4}\,.\end{split}$$

$$\begin{array}{r}\{\begin{array}{l}{\Vert \varphi \Vert }_{L^{\infty }({\Omega }_T)}\le C,\\ {\int }_0^{\infty }{\Vert \nabla \Delta \varphi \Vert }^{2}dt\le C(\lambda ),\\ {\Vert {\Delta }^2\varphi \Vert }_{L^{\infty }([0,\infty [;L^2(\Omega ))}\le C{\lambda }^{-C},\end{array}\end{array}$$

$$\begin{split}\phi_{0},\phi_{\epsilon}\in{\color[rgb]{0,0,0}C^{1}(0,T;W^{k,\infty}(\Omega^{\epsilon}))}\quad\text{for a }k\geq 4\,.\end{split}$$

$$\begin{array}{r}{\xi }_{\varphi }^k,{\xi }_w^k\in W^{1,\infty }(Y^1)\text{for all}1\le k\le d.\end{array}$$

$$\begin{array}{r}\{\begin{array}{ll}(T1)& T_{ϵ}u=u\text{a.e. in }{\Omega }^{ϵ},\\ (T2)& {\int }_{\Omega (ϵk_0)}{\left|T_{ϵ}u\right|}^{p}d\mathbf{x}\le k_1{\int }_{{\Omega }^{ϵ}}{\left|u\right|}^{p}d\mathbf{x},\\ (T3)& {\int }_{\Omega (ϵk_0)}{\left|D(T_{ϵ}u)\right|}^{p}d\mathbf{x}\le k_2{\int }_{{\Omega }^{ϵ}}{\left|Du\right|}^{p}d\mathbf{x},\end{array}\end{array}$$

$$\begin{split}\frac{\partial\mu_{0}(\phi_{0}({\bf x}))}{\partial x_{k}}=\begin{cases}0&\text{{\color[rgb]{0,0,0} appearing in the cell problem depending on $\Omega\times Y$}}\,,\\ \frac{\partial\mu_{0}(\phi_{0})}{\partial x_{k}}&\text{{\color[rgb]{0,0,0} appearing on the macroscale $\Omega$ (after averaging over $Y$)}}\,,\end{cases}\end{split}$$

$$\begin{array}{r}\{\begin{array}{ll}{\theta }_1\frac{\partial {\varphi }_0}{\partial t}=\mathrm{div}\left({\hat{\mathrm{M}}}_{\varphi }/\lambda \nabla f({\varphi }_0)\right)+\frac{\lambda }{{\theta }_1}\mathrm{div}\left({\hat{\mathrm{M}}}_w\nabla \left(\mathrm{div}\left(\hat{\mathrm{D}}\nabla {\varphi }_0\right)\right)\right)& \text{in }{\Omega }_T,\\ {\nabla }_n{\varphi }_0=\mathbf{n}\cdot \nabla {\varphi }_0=0& \text{on }\partial {\Omega }_T,\\ {\nabla }_n\Delta {\varphi }_0=0& \text{on }\partial {\Omega }_T,\\ {\varphi }_0(\mathbf{x},0)=\psi (\mathbf{x})& \text{in }\Omega ,\end{array}\end{array}$$

$$\begin{array}{r}\{\begin{array}{ll}{\mathrm{d}}_{ik}& :=\frac{1}{|Y|}{\sum }_{j=1}^d{\int }_{Y^1}\left({\delta }_{ik}-{\delta }_{ij}\frac{\partial {\xi }_{\varphi }^k}{\partial y_j}\right)d\mathbf{y},\\ {\mathrm{m}}_{ik}^{\varphi }& :=\frac{1}{|Y|}{\sum }_{j=1}^d{\int }_{Y^1}\left({\mathrm{m}}_{ik}-{\mathrm{m}}_{ij}\frac{\partial {\xi }_{\varphi }^k}{\partial y_j}\right)d\mathbf{y},\\ {\mathrm{m}}_{ik}^w(\mathbf{x})& :=\frac{1}{|Y|}{\sum }_{j=1}^d{\int }_{Y^1}\left({\mathrm{m}}_{ik}-{\mathrm{m}}_{ij}\frac{\partial {\xi }_w^k}{\partial y_j}\right)d\mathbf{y},\end{array}\end{array}$$

$$\begin{split}\xi_{w}^{k}:\,\,\begin{cases}-\sum_{i,j,k=1}^{d}\frac{\partial}{\partial y_{i}}\left({\rm m}_{ik}-{\rm m}_{ij}\frac{\partial\xi^{k}_{w}}{\partial y_{j}}\right)\\ \qquad\qquad=-\sum_{k,i,j=1}^{d}\frac{\partial}{\partial y_{i}}\left({\rm m}_{ik}-{\rm m}_{ij}\frac{\partial\xi^{k}_{\phi}}{\partial y_{j}}\right)&\textrm{ in }Y^{1}\,,\\ \sum_{i,j,k=1}^{d}{\rm n}_{i}\Bigl{(}\left({\rm m}_{ij}\frac{\partial\xi^{k}_{w}}{\partial y_{j}}-{\rm m}_{ik}\right)\\ \qquad\qquad+\left({\rm m}_{ik}-{\rm m}_{ij}\frac{\partial\xi^{k}_{\phi}}{\partial y_{j}}\right)\Bigr{)}=0&\textrm{ on }I_{Y}:=\partial Y^{1}\cap\partial Y^{2}\,,\\ \xi^{k}_{w}({\bf y})\textrm{ is $Y$-periodic and ${\mathcal{M}}_{Y^{1}}(\xi^{k}_{w})=0$,}\end{cases}\\ \xi_{\phi}^{k}:\,\,\begin{cases}-\sum_{i,j=1}^{d}\frac{\partial}{\partial y_{i}}\left(\delta_{ik}-\delta_{ij}\frac{\partial\xi^{k}_{\phi}}{\partial y_{j}}\right)=0&\textrm{ in }Y^{1}\,,\\ \sum_{i,j=1}^{d}{\rm n}_{i}\left(\delta_{ij}\frac{\partial\xi^{k}_{\phi}}{\partial y_{j}}-\delta_{ik}\right)=0&\textrm{ on }I_{Y}\,,\\ \xi^{k}_{\phi}({\bf y})\textrm{ is $Y$-periodic and ${\mathcal{M}}_{Y^{1}}(\xi^{k}_{\phi})=0$,}\end{cases}\end{split}$$

$$\begin{array}{rl}{\Vert E_{ϵ}^w(\cdot ,t)\Vert }_{L^2({\Omega }^{ϵ})}^2& +c(m,\lambda ,\kappa ){\int }_0^t{\Vert {\mathcal{A}}_{ϵ}{E}_{ϵ}^w(\cdot ,s)\Vert }_{L^2({\Omega }^{ϵ})}^{2}ds\\ & \le {ϵ}^{1/2}C(T,\Omega ,m,\kappa ,\lambda ),\\ {\Vert E_{ϵ}^{\varphi }(\cdot ,t)\Vert }_{H^1({\Omega }^{ϵ})}& \le {ϵ}^{1/4}C(T,\Omega ,m,\kappa ,\lambda ),\end{array}$$

$$\begin{split}\begin{cases}\quad{\mathcal{A}}_{0}&=-\sum_{i,j=1}^{d}\frac{\partial}{\partial y_{i}}\left(\delta_{ij}\frac{\partial}{\partial y_{j}}\right)\,,\\ \quad{\mathcal{A}}_{1}&=-\sum_{i,j=1}^{d}\biggl{[}\frac{\partial}{\partial x_{i}}\left(\delta_{ij}\frac{\partial}{\partial y_{j}}\right)+\frac{\partial}{\partial y_{i}}\left(\delta_{ij}\frac{\partial}{\partial x_{j}}\right)\biggr{]}\,,\\ \quad{\mathcal{A}}_{2}&=-\sum_{i,j=1}^{d}\frac{\partial}{\partial x_{j}}\left(\delta_{ij}\frac{\partial}{\partial x_{j}}\right)\,,\\ \quad{\mathcal{B}}_{0}&=-\sum_{i,j=1}^{d}\frac{\partial}{\partial y_{i}}\left({\rm m}_{ij}\frac{\partial}{\partial y_{j}}\right)\,,\\ \quad{\mathcal{B}}_{1}&=-\sum_{i,j=1}^{d}\biggl{[}\frac{\partial}{\partial x_{i}}\left({\rm m}_{ij}\frac{\partial}{\partial y_{j}}\right)+\frac{\partial}{\partial y_{i}}\left({\rm m}_{ij}\frac{\partial}{\partial x_{j}}\right)\biggr{]}\,,\\ \quad{\mathcal{B}}_{2}&=-\sum_{i,j=1}^{d}\frac{\partial}{\partial x_{j}}\left({\rm m}_{ij}\frac{\partial}{\partial x_{j}}\right)\,.\end{cases}\end{split}$$

$$\begin{array}{rl}{\zeta }^{ϵ}& \approx {\zeta }_0(\mathbf{x},\mathbf{y},t)+ϵ{\zeta }_1(\mathbf{x},\mathbf{y},t)+{ϵ}^2{\zeta }_2(\mathbf{x},\mathbf{y},t),\text{for}\zeta \in \left\{w,\varphi \right\},\end{array}$$

$$\begin{array}{r}f({\varphi }^{ϵ})\approx f({\varphi }_0)+f^{\prime }({\varphi }_0)({\varphi }^{ϵ}-{\varphi }_0)+\mathcal{O}\left(({\varphi }^{ϵ}-{\varphi }_0{)}^2\right),\end{array}$$

$$\begin{array}{r}\mathcal{O}({ϵ}^{-2}):\{\begin{array}{ll}{\mathcal{B}}_0\left[\lambda w_0+1/\lambda f({\varphi }_0)\right]=0& \text{in }{Y}^1,\\ \text{no flux b.c.},& \\ \text{w0 is Y1-periodic},& \\ {\mathcal{A}}_0{\varphi }_0=0& \text{in }{Y}^1,\\ {\nabla }_n{\varphi }_0=0& \text{on }\partial Y^1\cap \partial Y^2,\\ \text{ϕ0 is Y1-periodic},& \end{array}\end{array}$$

$$\begin{array}{r}\mathcal{O}({ϵ}^{-1}):\\ \{\begin{array}{ll}{\mathcal{B}}_0\left[\lambda w_1+1/\lambda f^{\prime }({\varphi }_0){\varphi }_1\right]=-{\mathcal{B}}_1\left[\lambda w_0+1/\lambda f({\varphi }_0)\right]& \text{in }{Y}^1,\\ \text{no flux b.c.},& \\ \text{w1 is Y1-periodic},& \\ {\mathcal{A}}_0{\varphi }_1=-{\mathcal{A}}_1{\varphi }_0& \text{in }{Y}^1,\\ {\nabla }_n{\varphi }_1=0& \text{on }\partial Y^1\cap \partial Y^2,\\ \text{ϕ1 is Y1-periodic},& \end{array}\end{array}$$

$$\begin{split}\mathcal{O}(\epsilon^{0}):\hskip 298.75394pt\\ \begin{cases}{\mathcal{B}}_{0}\left[\lambda w_{2}+\frac{1}{\lambda}\left(\frac{1}{2}f^{\prime\prime}(\phi_{0})\phi_{1}^{2}+f^{\prime}(\phi_{0})\phi_{2}\right)\right]&\\ \qquad=-\left({\mathcal{B}}_{2}\left[\lambda w_{0}+1/\lambda f(\phi_{0})\right]+{\mathcal{B}}_{1}\left[\lambda w_{1}{\color[rgb]{0,0,0}+}1/\lambda f^{\prime}(\phi_{0})\phi_{1}\right]\right)&\\ \qquad\qquad\qquad\quad-\partial_{t}\left(-\Delta\right)^{-1}w_{0}&\textrm{in }Y^{1}\,,\\ \quad\textrm{no flux b.c.}\,,\\ \quad\textrm{$w_{2}$ is $Y^{1}$-periodic}\,,\\ {\mathcal{A}}_{0}\phi_{2}=-{\mathcal{A}}_{2}\phi_{0}-{\mathcal{A}}_{1}\phi_{1}+w_{0}&\textrm{in }Y^{1}\,,\\ \quad\nabla_{n}\phi_{2}=g_{\epsilon}&\textrm{on }\partial Y^{1}\cap\partial Y^{2}\,,\\ \quad\textrm{$\phi_{2}$ is $Y^{1}$-periodic}\,.\end{cases}\end{split}$$

Full text

{AMS}

68Q25, 68R10, 68U05

1 Introduction

We consider the well-accepted Cahn-Hilliard/diffuse-interface formulation [9, 57] for studying the evolution of interfaces between different phases. Its broad applicability together with increasing computational power have enabled its use to new and increasingly complex scientific and engineering problems such as the computation of transport equations in porous media [42] which represents a numerically very demanding, high-dimensional multiscale problem [26]. The purpose of the present work is to rigorously and systematically provide a reliable effective macroscopic description of how multiple phases invade strongly heterogeneous media, such as porous materials for instance.

The cornerstone of phase-field models is the abstract energy density

[TABLE]

where $\phi:=\frac{c_{\beta}}{c_{\alpha}+c_{\beta}}$ is a reduced order parameter representing the fraction of species of type $\beta$ in a binary solution containing species $\alpha$ and $\beta$ with number densities $c_{\alpha}$ and $c_{\beta}$, respectively. The gradient term {\color[rgb]{0,0,0}\lambda/2}\left|\nabla\phi\right|^{2} penalizes the interfacial area between these phases, and $F$ is defined as the general (Helmholtz) free energy density $F(\phi):=U-TS$ , where $U$ is the internal energy, $T$ is the temperature and $S$ is the entropy. The parameter $\lambda>0$ is proportional to the interfacial width and leads to the appearance of smooth interface.

Important examples of this formulation include the regular solution theory which has been applied successfully in a wide spectrum of scientific and technological contexts such as ionic melts [23], water sorption in porous solids [5], and micellization in binary surfactant mixtures [25]. The key quantity in this theory is the so-called regular solution energy density (also known as the Flory-Huggins energy density [18]) $F(\phi):=R(\phi)-TS_{I}(\phi)$ , where $S_{I}(\phi):=-k_{B}\left[\phi{\rm ln}\,\phi-(1-\phi){\rm ln}\,(1-\phi)\right]$ is the entropy of mixing for ideal solutions and the regular solution term $R(\phi):=z\omega\phi(1-\phi)$ accounts for the interaction energy between different species. The variable $z$ is the coordination number defining the number of bonds of $\beta$ with neighbouring species. $\omega:=\epsilon_{\alpha\alpha}+\epsilon_{\beta\beta}-2\epsilon_{\alpha\beta}$ is the interaction energy parameter accounting for the minima $\epsilon_{\alpha\alpha}$, $\epsilon_{\beta\beta}$, and $\epsilon_{\alpha\beta}$ of interaction potentials which define attractive and repulsive forces between the species $\alpha$ and $\beta$.

Wetting phenomena, often studied using classical sharp-interface approximations, e.g. [43, 45, 44, 58], also enjoy a wide-spread use of phase-field modeling [41, 60, 59, 52, 53] even in the presence of complexities such as an electric field (so called electrowetting, e.g. [14, 34]). The reason for this is that classical sharp-interface models consider the fluid-fluid interface to be a sharp surface of zero thickness where quantities such as the fluid density are, in general, discontinuous, which leads to singularity formation for interfacial problems with topological transitions, e.g. the notorious contact line singularity [22] often cured with phenomenological approaches such as slip models. The phase-field/diffuse-interface approach relaxes the assumption of a sharp interface in line with the physics of the problem and in agreement with developments and applications in the field of statistical mechanics of liquids and in molecular simulations, with quantities varying smoothly but rapidly, and considers the interface to have a non-zero thickness, thus allowing a “natural” regularisation for singularities in interfacial problems with topological transitions.

Other applications include transport in electrochemical systems e.g. consisting of an electrolyte and an electrode [20], or immiscible flows [28, 37] under a polynomial free energy in the form of the classical double-well potential, i.e., $W(\phi):=\frac{1}{4}(1-\phi^{2})^{2}$ are relevant applications. Phase-field energy functionals are also of interest in image processing such as inpainting, see e.g. [7].

Our formal derivation of upscaled phase-field equations is valid for general free energies but the subsequent rigorous derivation of error estimates is based on free energies of the following form.

Polynomial Class (PC): Admissible free energy densities $F$ in (1) are polynomials of order {\color[rgb]{0,0,0}2r}, i.e.,

[TABLE]

with $f(u)=F^{\prime}(u)$ vanishing at $u=0$, that is,

[TABLE]

*where the leading coefficient of both $F$ and $f$ is positive, i.e., $a_{2r-1}=2rb_{2r}>0$. *

Temam [56] established well-posedness of the Cahn-Hilliard equation for free energies of class (PC). In computations, one often replaces the regular solution energy density, composed of $R$ and $S_{I}$ defined above, by the polynomial double-well potential $W(\phi)$.

In difference to [51], we provide here an upscaling strategy that is valid for general homogeneous free energy densities by making use of a Taylor expansion of the free energy density at the effective upscaled solution. This serves also as a general methodology for the homogenization of nonlinear problems. Moreover, to the best of our knowledge, we present here for the first time, error estimates between the solution of the microscopic phase-field equations solved in a periodic porous medium and the solution of the correspondingly homogenized/upscaled equations by Theorem 1 below. In the remaining part of this section, we introduce the basic equations describing interfacial dynamics in a homogeneous environment and subsequently in a periodic porous medium.

(a) Homogeneous domains $\Omega$. In the Ginzburg-Landau/Cahn-Hilliard formulation, the total energy is defined by $E(\phi):=\int_{\Omega}e(\phi)\,d{\bf x}$ with density (1) on a bounded domain $\Omega\subset\mathbb{R}^{d}$ with smooth boundary $\partial\Omega$ and $1\leq d\leq 3$ denotes the spatial dimension. It is well accepted that thermodynamic equilibrium can be achieved by minimizing the energy $E$, frequently supplemented by a wetting boundary contribution $\int_{\partial\Omega}g({\bf x})\,d{\bf x}$ for $g({\bf x})\in H^{3/2}(\partial\Omega)$. The wetting property of pore walls can be characterised by

[TABLE]

where $C_{h}$ is the Cahn number $\frac{\lambda}{L}$, $L$ the macroscopic length scale, and $\gamma=\frac{2\sqrt{2}\phi_{e}}{3\sigma_{lg}}$, $\sigma_{lg}$ is the liquid-gas surface tension, and $\phi_{e}$ is the local equilibrium limiting value of $F$, see [46]. For simplicity, we set subsequently $g=0$ and hence assume walls with neutral wetting characteristics, i.e., walls inducing a contact angle of 90 degrees. A widely used minimization over time forms the $H^{-1}$-gradient flow with respect to $E(\phi)$, i.e.,

[TABLE]

where $\Omega_{T}:=\Omega\times]0,T[$, $\partial\Omega_{T}:=\partial\Omega\times]0,T[$, $\phi$ satisfies the initial condition $\phi({\bf x},0)=\psi({\bf x})$, and $\hat{\rm M}=\left\{{\rm m}_{ij}\right\}_{1\leq i,j\leq d}$ denotes a symmetric and positive definite mobility tensor. Throughout the article we write $]a,b[$ for open intervals with $a,b\in\mathbb{R}$ and $a<b$. The gradient flow (4) is weighted by the mobility tensor $\hat{\rm M}$, and is referred to as the Cahn-Hilliard equation. This equation is a model prototype for interfacial dynamics, e.g. [17], and phase transformation, e.g. [9], under homogeneous Neumann boundary conditions, i.e., $g=0$, and free energy densities $F$.

We recall that the integrated energy density (1) dissipates along solutions of (4), that means, $E(\phi(\cdot,t))\leq E(\phi(\cdot,0))=:E_{0}$. This follows immediately after differentiating $E(\phi)$ with respect to time and using (4) for $g=0$.

There is also an interesting connection between the Cahn-Hilliard/phase-field equation and the free-boundary value problem known as the Mullins-Sekerka problem [32] or the two-phase Hele-Shaw problem [21]. The Hele-Shaw problem plays a crucial role for deriving more regular solutions of the Cahn-Hilliard equation (4), see [2]. Inspired by the formal derivation by Pego [40], it was rigorously verified later on in [2, 54] that the chemical potential

[TABLE]

satisfies for an evolving interfacial front $\Gamma_{t}$ with initial condition $\Gamma_{00}$ in the limit $\lambda\to 0$ for $t\in[0,T]$ the following

[TABLE]

where $\sigma=\int_{-1}^{1}\left(\frac{1}{2}\int_{0}^{s}f(r)\,dr\right)^{1/2}\,ds$ is the interfacial tension, $\kappa$ the mean curvature, $v$ the normal velocity of the interface $\Gamma_{t}$, ${\bf n}$ the unit outward normal to either $\partial\Omega$ or $\Gamma_{t}$, and $\left[{\bf n}\cdot\nabla\mu\right]_{\Gamma_{t}}:={\bf n}\cdot\nabla\mu^{+}-{\bf n}\cdot\nabla\mu^{-}$ where \mu^{+}:=\mu\,\bigr{|}_{\Omega^{+}_{t}} and \mu^{-}:=\mu\,\bigl{|}_{\Omega^{-}_{t}} and $\Omega^{+}_{t}$ and $\Omega^{-}_{t}$ denote the exterior and interior of $\Gamma_{t}$ in $\Omega$. Herewith, we also have $\phi\to\pm 1$ in $\Omega^{\pm}_{t}$ for all $t\in[0,T]$ as $\lambda\to 0$. Finally, the derivation of convergence rates (Theorem 1 below) requires higher regularity of solutions of the Cahn-Hilliard equation (Assumption C below) than available in [2, 16], which require the existence of global in time solutions of the sharp interface limit (6).

(b) Heterogeneous/perforated domains $\Omega^{\epsilon}$. Our main study concentrates on (1) in perforated domains $\Omega^{\epsilon}\subset\mathbb{R}^{d}$ instead of a homogeneous domain $\Omega\subset\mathbb{R}^{d}$. The dimensionless variable $\epsilon>0$ defines the heterogeneity $\epsilon=\frac{\ell}{L}$ where $\ell$ represents the characteristic pore size and $L$ is the macroscopic length of the porous medium, see Figure 1. Hence, the porous medium is defined by a reference pore/cell $Y:=[0,\ell_{1}]\times[0,\ell_{2}]\times\dots\times[0,\ell_{d}]$. For simplicity, we set $\ell_{1}=\ell_{2}=\dots=\ell_{d}=1$. The pore and the solid phase of the medium are denoted by $\Omega^{\epsilon}$ and $B^{\epsilon}$, respectively. These sets are defined by,

[TABLE]

where the subsets $Y^{1},\,Y^{2}\subset Y$ are such that $\Omega^{\epsilon}$ is a connected set. More precisely, $Y^{1}$ stands for the pore phase (e.g. liquid or gas phase in wetting problems), see Figure 1. Additionally, we define the macroscopic pore walls by $I_{\Omega}^{\epsilon}:=\partial\Omega^{\epsilon}\cap\partial B^{\epsilon}$ and the microscopic pore walls by $I_{Y}:=\partial Y^{1}\cap\partial Y^{2}$. Herewith, we can reformulate (4) for $g=0$ by the following microscopic porous media problem

[TABLE]

Our main objective is the derivation of error estimates for the difference between the upscaled/homogenized solution $\phi_{0}$ of (17) and the microscopic solution $\phi_{\epsilon}$ of (8) in order to have a qualitative and quantitative measure for the validity of the homogenized phase field formulation (17) (Theorem 1) obtained by passing to the limit $\epsilon\to 0$ in (8). This result will also provide a rigorous basis for the formal upscaling in [50]. The homogenized equation stated in Theorem 1 below allows for new analytical considerations such as a sharp interface study of the novel upscaled equation or establishing more regular solutions of Cahn-Hilliard/phase field equations as well as for new avenues in modelling. It ultimately leads to convenient, low-dimensional computational schemes which can be solved by well-known numerical methods developed for homogeneous domains.

In Section 2, we present basic notations and mathematical assumptions. The main results are summarized in Section 3 and subsequently justified in Sections 4 and 5. Conclusions and suggestions for further work are given in Section 6.

2 Mathematical preliminaries and notation

We recall the splitting formulation of the Cahn-Hilliard equation from [51] which builds the basis for our subsequent homogenization analysis. To this end, we set

[TABLE]

and identify $\phi=(-\Delta)^{-1}w$ in the $H^{2}_{E}(\Omega)$-sense, this means, we have for all $\varphi\in H^{2}_{E}(\Omega)$ that

[TABLE]

where $(\cdot,\cdot)$ denotes the standard $L^{2}$-scalar product. Herewith we can rewrite (8) (for simplicity stated here for $\Omega\subset\mathbb{R}^{d}$ instead of $\Omega^{\epsilon}$) for all $\varphi\in H^{2}_{E}(\Omega)$ as

[TABLE]

which reads in the classical sense

[TABLE]

In [35], the existence of a local solution $\phi\in H^{2}_{E}(\Omega)$ to equation (8) has been verified for $f\in C^{2}_{Lip}(\mathbb{R})$ and hence also to (9). Furthermore, Novick-Cohen [35] states necessary conditions for global existence while a proof based on Galerkin approximations and a priori estimates can be found in [56, Theorem 4.2, p. 155]. The advantage of (9) is that it allows to base our upscaling approach on well-known results from elliptic/parabolic homogenization theory [6, 11, 24, 30, 39, 62]. The splitting (9) slightly differs from the strategy of substituting the chemical potential, which is often applied for computational purposes, see [4], and which seems also more appropriate for other homogenization strategies such as periodic unfolding [12] or two-scale convergence [3, 33] for instance.

Next, we briefly summarize what, to the best of our knowledge, we believe to be the best available regularity results (Lemma 1 below) for the Cahn-Hilliard equation [2, 16]. These results depend on two assumptions:

Assumption A: **

• (A1)

$F\in C^{4}(\mathbb{R})$* satisfies $F(\pm 1)=0$ and $F>0$ elsewhere.*

• (A2)

$f(u)=F^{\prime}(u)$* satisfies for some finite $\alpha>2$ and positive constants $k_{i}>0$, $i=0,\dots,3$,*

[TABLE]

• (A3)

There exist constants $0<a_{1}\leq 1$, $a_{2}>0$, $a_{3}>0$ and $a_{4}>0$ such that for $b\in\mathbb{R}$

[TABLE]

It is straightforward to check that the classical double-well potential $F(x)=(x^{2}-1)^{2}/4$ satisfies Assumption A. The following characterization of the initial condition $\psi$ is also required for more regular solutions as derived in [2, 16] and stated in Lemma 1 below. We will frequently write $\left\|u\right\|$ for the $L^{2}$-norm of a function $u$.

Assumption B: There exist uniform constants $m_{0},\,\sigma_{j}>0$, $j=1,2,3$ such that

• (B1)

$-1<m_{0}:=\frac{1}{\left|\Omega\right|}\int_{\Omega}\psi({\bf x})\,d{\bf x}<1\,,$**

• (B2)

${\cal E}_{\lambda}(\psi):=\frac{\lambda}{2}\left\|\nabla\psi\right\|^{2}+\frac{1}{\lambda}\left\|F(\psi)\right\|_{L^{1}}\leq C\lambda^{-2\sigma_{1}}\,,$**

• (B3)

$\left\|\omega^{\lambda}\right\|_{H^{l}}:=\left\|-\lambda\Delta\psi+\frac{1}{\lambda}F(\psi)\right\|_{H^{l}}\leq C\lambda^{-\sigma_{2+l}}\,,\quad l=0,1\,,$**

*where $|\Omega|$ is the Lebesgue measure of $\Omega|$. *

Herewith, the following regularity result has been derived for homogeneous domains $\Omega$ in [2, 16].

Lemma 1**.**

(Regularity)* Let $f$ and $\psi$ satisfy the Assumption A and B, respectively. Moreover, we suppose that the Hele-Shaw/Mullins-Sekerka problem (6) has a global in time classical solution. Then, the solution $\phi$ of the Cahn-Hilliard equation (4) satisfies the estimates*

[TABLE]

*for all $\lambda\in]0,\kappa[$ and a family of smooth initial data $\left\{\psi^{\lambda}\right\}_{0<\lambda\leq 1}$ where $\kappa$ and $C$ are constants. Estimate (12)3 holds for $C>0$ large enough, if $\lim_{s\to 0^{+}}\left\|\nabla\partial_{t}\phi(s)\right\|\leq C\lambda^{-\kappa}$. *

Remark 1**.**

(Hele-Shaw)* Existence and uniqueness of classical solutions for the so-called single phase Hele-Shaw problem in bounded domains in $\mathbb{R}^{d}$ can be found for instance in [15, 31]. $\diamond$ *

We refer the interested reader to Refs. [2, 16] for a proof. Since we need slightly stronger regularity results than stated in Lemma 1 for the proof of error estimates (Theorem 1), we introduce the following well-accepted (see for instance [11])

Assumption C: * For smooth data, i.e., $\phi_{0}({\bf x},0),\phi_{\epsilon}({\bf x},0)\in C^{\infty}(\Omega^{\epsilon})$, $f\in C^{\infty}(\mathbb{R})$, and for $\Omega^{\epsilon}$ with Lipschitz boundary $\partial\Omega^{\epsilon}$ and hence the interface $I^{\epsilon}_{\Omega}$ is Lipschitz too, then the solutions $\phi_{0}$ of equation (17) and $\phi_{\epsilon}$ the solution of the microscopic equation (8) satisfy*

[TABLE]

Moreover, the correctors $\xi_{\phi}^{k}$ and $\xi_{w}^{k}$, which solve the cell problems (19), satisfy

[TABLE]

Let $T_{\epsilon}$ denote the extension operator, which extends the solutions $\phi^{\epsilon}$ and $w^{\epsilon}$ of (9) defined on the perforated domain $\Omega^{\epsilon}$ to the homogeneous domain $\Omega$. For convenience we denote these extensions by $\phi^{\epsilon}$ and $w^{\epsilon}$ and skip the extension operator $T_{\epsilon}$ most of the time. The existence of such an operator $T_{\epsilon}\,:\,W^{1,p}(\Omega^{\epsilon})\to W^{1,p}_{loc}(\Omega)$ for $\epsilon>0$ was established in [1] and $T_{\epsilon}$ is characterized by the following properties:

[TABLE]

for constants $k_{0},\,k_{1},\,k_{2}>0$. Hence, $T_{\epsilon}$ extends solutions defined on the pore space $\Omega^{\epsilon}$ to the whole domain $\Omega$.

3 Main results

Our main result, i.e., the upscaling/homogenization of general phase field equations (including the Cahn-Hilliard equation), is based on the following local property of the chemical potential.

Definition 1**.**

(Local Thermodynamic Equilibrium)* Let $\mu(\phi)=-\lambda\Delta\phi+\frac{1}{\lambda}f(\phi)$ be the chemical potential associated to the phase field free energy density (1). We say that the upscaled chemical potential
$\mu_{0}(\phi_{0})=-\lambda{\rm div}(\hat{\rm D}\nabla\phi_{0})+\frac{1}{\lambda}f(\phi_{0})=\lambda w_{0}+\frac{1}{\lambda}f(\phi_{0})$ is in local thermodynamic equilibrium (LTE) if and only if*

[TABLE]

*where $\phi_{0}(x)$ is the upscaled/slow variable, which is independent of the microscale ${\bf y}\in Y$ and which solves the upscaled phase field equation (17) below. *

Remark 2**.**

*Definition 1 systematically accounts for the problem specific slow (macroscopic) scale ${\bf x}\in\Omega$ and the fast (microscopic) scale ${\bf y}\in Y$. Intuitively, Definition 1 expresses the fact that the macroscopic variables are varying so slowly that their variations are not visible on the microscale. The local thermodynamic equilibrium characterization (16) is well accepted and appears in a wide range of applications, e.g. [8, 13, 29, 27]. $\diamond$ *

Definition 1 naturally appears in the upscaling of nonlinear problems and enables two essential features: a) The upscaled equations are of the same form as the microscopic formulation; b) (16) guarantees the well-posedness of arising cell problems which define effective transport coefficients. Recent examples in the context of ionic transport equations are [46, 47, 49, 48]. These considerations allow us to recall the following upscaling result from [50].

Upscaling Result (UR): (Effective macroscopic phase field equations) * Suppose that $\psi({\bf x})\in H^{2}_{E}(\Omega)$. For chemical potentials $\mu:=\nabla_{\phi}E(\phi)$, where $\nabla_{\phi}$ denotes the Fréchet derivative, being in local thermodynamic equilibrium as characterized by Definition 1, the microscopic porous media formulation (8) can be effectively approximated by the following macroscopic problem,*

[TABLE]

where $\theta_{1}:=\frac{\left|Y^{1}\right|}{\left|Y\right|}$ is the porosity and the porous media correction tensors $\hat{\rm D}:=\left\{{\rm d}_{ik}\right\}_{1\leq i,k\leq d}$, $\hat{\rm M}_{\phi}=\left\{{\rm m}^{\phi}_{ik}\right\}_{1\leq i,k\leq d}$ and \hat{\rm M}_{w}=\left\{{\color[rgb]{0,0,0}{\rm m}^{w}_{ik}}\right\}_{1\leq i,k\leq d} are defined by

[TABLE]

where $m_{ij}$ are elements of the mobility tensor. The corrector functions $\xi^{k}_{\phi}\in H^{1}_{per}(Y^{1})$ and $\xi^{k}_{w}\in L^{2}(\Omega;H^{1}_{per}(Y^{1}))$ for $1\leq k\leq d$ solve in the distributional sense the following reference cell problems

[TABLE]

*where $\delta_{ij}$ is the Kronecker delta function and ${\rm n}_{i}$ denotes the $i$-th component of the outward normal vector $\bf n$. *

Remark 3**.**

*i) For an isotropic mobility, i.e., $\hat{\rm M}:=m\hat{\rm I}$ where $\hat{\rm I}$ is the identity matrix, it follows that $\xi^{k}_{w}=\xi^{k}_{\phi}$, and hence both $\xi_{w}^{k}$ and $\xi_{\phi}^{k}$ solve classical elliptic cell problems, see Figure 2.

ii) The local thermodynamic equilibrium property (Definition 1) of the macroscopic chemical potential $\mu_{0}$ enables the derivation of the well-posed cell problem (19)1 for $\xi_{w}^{k}$. $\diamond$ *

The next result characterizes qualitatively the homogenized phase field equations (17) with the help of error estimates.

Theorem 1**.**

(Error estimates)* Let $\phi^{\epsilon}$ be a solution of (8), or equivalently $\phi^{\epsilon}$ and $w^{\epsilon}$ solve the splitting formulation (9). Suppose that Assumption C holds. Moreover, the domain boundaries $\partial\Omega^{\epsilon}$ and interfaces $I^{\epsilon}_{\Omega}:=\partial\Omega^{\epsilon}\cap\partial B^{\epsilon}$ shall be Lipschitz111see [10] for instance. Let $\hat{\rm M}=m\hat{\rm I}$ be an isotropic mobility with $\hat{\rm I}$ representing the identity matrix. If the free energy $F$ is polynomial of class (PC), then the error variables ${\rm E}^{\phi}_{\epsilon}:=\phi^{\epsilon}-(\phi_{0}+\epsilon\phi_{1})\,,{\rm E}^{w}_{\epsilon}:=w^{\epsilon}-(w_{0}+\epsilon w_{1})$, where $w_{1}:=-\sum_{k=1}^{d}\xi^{k}_{w}({\bf y})\frac{\partial w_{0}}{\partial x_{k}}({\bf x},t)$ and $\phi_{1}:=-\sum_{k=1}^{d}\xi^{k}_{\phi}({\bf y})\frac{\partial\phi_{0}}{\partial x_{k}}({\bf x},t)$, satisfy for $0\leq t\leq T$ and $0<T<\infty$ the following estimates*

[TABLE]

*where $c(m,\lambda,\kappa)$ and $C(T,\Omega,m,\kappa,\lambda)$ are constants independent of $\epsilon$. *

Remark 4**.**

*We note that the proof of the above Theorem 1 does not take the behaviour in the boundary region into account by solely applying a smooth enough truncation. This leads for linear elliptic equations to the by now classical convergence rate $\epsilon^{1/2}$, e.g. [10, 61]. However, in recent attempts [55, 38], the authors can improve the convergence rates with the help of operator estimates with a resulting rate $\epsilon$. We note, that our estimates in (20) are derived based on the classical method but due to the fourth order operator, we end up with the slightly lower rate $\epsilon^{1/4}$, albeit under the generally required strong regularity Assumption C. The strongest regularity result currently available seems to be the estimates stated in Lemma 1. *

To the best of our knowledge, this is the first error quantification in terms of convergence rates with respect to the heterogeneity $\epsilon$ of the porous media approximation (17) for phase field equations. The estimates (20) imply convergence of solutions $\phi^{\epsilon}$ of the microscopic formulation (8) to solutions $\phi_{0}$ of the upscaled problem (17) for a vanishing heterogeneity parameter based on the regularity Assumption C.

4 Formal derivation of upscaled equations

For convenience, we first recall here the formal derivation of the effective macroscopic phase field equation from [50]. For the micro-scale variable $\frac{{\bf x}}{\epsilon}=:{\bf y}\in Y$ it holds that, \frac{\partial f_{\epsilon}({\bf x})}{\partial x_{i}}=\frac{1}{\epsilon}\frac{\partial f}{\partial y_{i}}({\bf x},{\bf x}/\epsilon)+\frac{\partial f}{\partial x_{i}}({\bf x},{\bf x}/\epsilon)\,,\textrm{ and }\nabla f_{\epsilon}({\bf x})={\color[rgb]{0,0,0}\frac{1}{\epsilon}\nabla_{y}f({\bf x},{\bf x}/\epsilon)+\nabla_{x}f({\bf x},{\bf x}/\epsilon)}\,,  where $f_{\epsilon}({\bf x})=f({\bf x},{\bf y})$ is an arbitrary function depending on two variables ${\bf x}\in\Omega$, ${\bf y}\in Y$. Hence, we have

[TABLE]

Herewith, we can define ${\mathcal{A}}_{\epsilon}:=\epsilon^{-2}{\mathcal{A}}_{0}+\epsilon^{-1}{\mathcal{A}}_{1}+{\mathcal{A}}_{2}$ and analogously ${\mathcal{B}}_{\epsilon}$. Hence, it holds for the Laplace operator that {\color[rgb]{0,0,0}-\Delta f_{\epsilon}({\bf x})={\mathcal{A}}_{\epsilon}f({\bf x},{\bf y})}. In order to deal with the multiscale nature of strongly heterogeneous environments [24, 46, 49, 47], the following, formal asymptotic expansions are used,

[TABLE]

where higher order terms are neglected. Before we can insert (22) into the microscopic formulation (9), we need to approximate the derivative of the nonlinear homogeneous free energy $f:=F^{\prime}$ by a Taylor expansion of the form

[TABLE]

where $\phi_{0}$ stands for the leading order term in (22).222 Here, we allow for general free energy densities in difference to the subsequent rigorous derivation of error estimates which is based on energy densities of the polynomial class (PC).

Using (22) and (23) in (9) with $\nabla_{n}\phi=g$ and with (21), we get the following sequence of problems,

[TABLE]

[TABLE]

[TABLE]

Problem (24) immediately suggests, based on classical homogenization theory [6], that $\phi_{0}$ is independent of the microscale ${\bf y}$. This and the linear structure of (25) allow for the following ansatz for $w_{1}$ and $\phi_{1}$, i.e.,

[TABLE]

Plugging (27) into (25)2 gives an equation for $\xi^{k}_{w}$ and $\xi^{k}_{\phi}$. The resulting equation for $\xi^{k}_{\phi}$ can be immediately written for $1\leq k\leq d$ as,

[TABLE]

To study (25)1, we first rewrite ${\mathcal{B}}_{0}\left[f^{\prime}(\phi_{0})\phi_{1}\right]$ and ${\mathcal{B}}_{1}f(\phi_{0})$ as follows

[TABLE]

Rewriting $w_{1}$ and $w_{0}$ in the same way and using (27) leads then to

[TABLE]

in $Y^{1}$. Next, due to local thermodynamic equilibrium property of the upscaled chemical potential $\mu_{0}(\phi_{0})$ as defined in Definition 1, we have on the level of the reference cell $Y$,

[TABLE]

Entering with (31) into (30) finally gives the reference cell problem for $\xi^{k}_{w}$, $1\leq k\leq d$ for given $\xi^{k}_{\phi}$

[TABLE]

Finally, we consider the last problem (26). Standard existence and uniqueness results (Fredholm alternative/Lax-Milgram) guarantee solvability after validating that the right hand side in (26) is zero as an integral over $Y^{1}$. This means, $-\sum_{i,j=1}^{d}\int_{Y^{1}}\frac{\partial}{\partial x_{i}}\left(\delta_{ij}\left(\frac{\partial\phi_{1}}{\partial y_{j}}+\frac{\partial\phi_{0}}{\partial x_{j}}\right)\right)\,d{\bf y}-\tilde{g}_{0}=\left|Y^{1}\right|w_{0}$ , where

[TABLE]

which leads to

[TABLE]

The inhomogeneous Neumann boundary condition $\tilde{g}_{0}$ accounts for pore walls $\partial Y^{1}_{w_{1}}\cup\partial Y^{1}_{w_{1}}=\partial Y^{1}$ showing two specific wetting properties characterised by the parameters $a_{1}$ and $a_{2}$ specifying the walls $\partial Y^{1}_{w_{1}}$ and $\partial Y^{1}_{w_{2}}$, respectively. The upscaling result (17) is stated for neutral wetting characteristics of the pore walls, i.e., $\tilde{g}=0$.

(33) suggests to define a porous media correction tensor $\hat{\rm D}:=\left\{{\rm d}_{ik}\right\}_{1\leq i,k\leq d}$ by

[TABLE]

Equations (33) and (34) represent the upscaled equation for $\phi_{0}$, i.e., $-\Delta_{\hat{\rm D}}\phi_{0}:=-{\rm div}\left(\hat{\rm D}\nabla\phi_{0}\right)=\theta_{1}w_{0}+\frac{1}{|Y|}\tilde{g}_{0}$.

The upscaled equation for $w$ is again a result of the Fredholm alternative, i.e., a solvability criterion on equation (26)1. We require,

[TABLE]

The first two terms in (35) can be rewritten by,

[TABLE]

where the effective tensor $\hat{\rm M}_{w}=\left\{{\rm m}^{w}_{ik}\right\}_{1\leq i,k\leq d}$ is defined by

[TABLE]

The third term in (35) becomes

[TABLE]

where the last term in (37) disappears after integrating by parts. The first term on the right-hand side of (37) can be rewritten with the help of the chain rule $\frac{\partial^{2}f(\phi_{0})}{\partial x_{k}\partial x_{j}}=f^{\prime\prime}(\phi_{0})\frac{\partial\phi_{0}}{\partial x_{k}}\frac{\partial\phi_{0}}{\partial x_{j}}+f^{\prime}(\phi_{0})\frac{\partial^{2}\phi_{0}}{\partial x_{k}\partial x_{j}}$ , as follows

[TABLE]

to which we add the term $-{\cal B}_{2}f(\phi_{0})$. Herewith, we define a tensor $\hat{\rm M}_{\phi}=\left\{{\rm m}^{\phi}_{ij}\right\}_{1\leq i,k\leq d}$ by

[TABLE]

which allows us to write

[TABLE]

These considerations together with the identity

[TABLE]

where the last term subsequently disappears due to $Y$-periodicity, finally lead after integration over the microscale $Y$ to the following effective equation for $\phi_{0}$, i.e.,

[TABLE]

The solvability of (39) follows along with the arguments in [35] via a local Lipschitz argument, or via a Galerkin approximation and a priori estimates as developed in [56, Theorem 4.2, p. 155].

5 Proof of Theorem 1

For the derivation of the error estimates (20), we work with the splitting formulation introduced in [51] and summarized in (9). We extend the derivation of error estimates for second order problems [11, Theorem 6.3, Section 6.2 & Section 7.2] and [46] to the fourth order phase field problems studied here. Hence, we compare the solution of the microscopic porous media formulation (8) with the solution of the effective upscaled porous media formulation (17). As in [11, 46], we introduce the error variables

[TABLE]

here extended by the error function $E^{f}_{\epsilon}$. The first goal is to determine the variables ${\rm F}_{\epsilon}^{\iota}$ and ${\rm G}_{\epsilon}^{\iota}$ which allow us to write the equations for the errors ${\rm E}_{\epsilon}^{\iota}$ for $\iota\in\left\{w,\phi\right\}$ as follows

[TABLE]

With the definitions (21) we can rewrite the first term on the right-hand side in (41)1 and the term on the left-hand side in (41)3 as follows

[TABLE]

The relations (21) and (42) together with the sequence of problems (24), (25), and (26) define the terms ${\rm R}^{\phi}_{\epsilon}$, ${\rm F}^{w}_{\epsilon}$ and ${\rm F}^{\phi}_{\epsilon}$ by

[TABLE]

since ${\cal B}_{\epsilon}=-m\Delta$. The definitions in (43) are a consequence of the two identities

[TABLE]

where the terms in the braces vanish due to the microscopic and homogenized equations and these terms represent the first and the second term in the error equation (41)1 and (41)3.

The inhomogeneities in the boundary conditions in (41) satisfy

[TABLE]

since the boundary conditions imply $\nabla_{n}(\iota^{\epsilon}-\iota_{0})=0$ for $\iota\in\left\{w,\phi\right\}$. Under the given regularity (Assumption C) of boundaries and interfaces, i.e., $\partial\Omega\in C^{\infty}$ as well as $I_{Y}\in C^{\infty}$, we obtain that $\xi_{\phi}^{k},\,\xi_{w}^{k}\in W^{1,\infty}$ by classical regularity theory for elliptic problems [19]. We note that the above regularity requirements on $\partial\Omega$ and $I_{Y}$ are not necessarily sharp but this question is beyond the scope of this work. Hence, elliptic theory allows us also to estimate (41)3 by

[TABLE]

where we subsequently justify the uniform boundedness of ${\rm F}^{\phi}_{\epsilon}$ in $L^{2}(\Omega^{\epsilon})$ and the $\epsilon$-dependent bound on ${\rm G}_{\epsilon}^{\phi}\in H^{1/2}(\partial\Omega^{\epsilon})$.

Next, we derive bounds for the terms on the right-hand side, i.e., (43). Thanks to Assumption C, it holds that

[TABLE]

for a constant $C>0$ independent of $\epsilon$. Analogously, we obtain the bound

[TABLE]

We show the basic steps to bound the second last term in (48). To this end, we first note that ${\cal B}_{2}[f^{\prime}(\phi_{0})\phi_{1}]=-m\Delta_{x}[f^{\prime}(\phi_{0})(\boldsymbol{\xi}_{\phi}\cdot\nabla_{x})\phi_{0}]$ from which we can identify the three most challenging terms to estimate:

[TABLE]

where ${\bf\xi}_{\phi}$ is the vector consisting of the corrector elements $\xi_{\phi}^{k}$ defined in (28). In order to bound the terms containing $\phi_{0}$ we use Assumption C. The factor $f^{\prime}(\phi_{0})$ and associated derivatives in (49) can be bounded by the fact that $f(s)$ is a polynomial of order $2p-1$ with characterization (PC), i.e.,

[TABLE]

We also note that $-\Delta_{x}\xi^{k}_{\phi}=-1/\epsilon^{2}\Delta_{y}\xi^{k}_{\phi}=0$ in view of the cell/corrector problem for $\xi^{k}_{\phi}$, $1\leq k\leq d$.

Finally, the remaining term ${\rm R}_{\epsilon}^{\phi}$ first decomposes in view of Assumption C as follows,

[TABLE]

and hence

[TABLE]

where $\left\|\Delta{\rm E}^{\phi}_{\epsilon}\right\|\leq\left\|{\rm E}^{w}_{\epsilon}\right\|+\epsilon\left\|{\rm F}^{\phi}_{\epsilon}\right\|$ thanks to (41)3

and $I_{\phi}:=]\underline{\phi},\overline{\phi}[$ is the interval defined by the smallest real root $\underline{\phi}$ and $\overline{\phi}$ the largest real root of the polynomial free energy $F$ characterized by (2).

In order to control the boundary contributions (45), we apply a standard argument [6, 11] based on a cut-off function $\chi^{\epsilon}$ which is defined as follows,

[TABLE]

We first look at ${\rm G}_{\epsilon}^{\phi}\,.$

For $\eta^{\phi}_{\epsilon}:=\chi^{\epsilon}{\rm G}^{\phi}_{\epsilon}$, we show that $\eta^{\phi}_{\epsilon}\in H^{1}(\Omega^{\epsilon})$ and

[TABLE]

where $U^{\epsilon}$ is the support of $\eta^{\phi}_{\epsilon}$ and forms a neighbourhood of $\partial\Omega^{\epsilon}$ of thickness $2\epsilon$. The regularity properties of $\xi^{k}_{w}$, $\xi_{\phi}^{k}$ and $\chi^{\epsilon}$ allow us to control $\eta_{\epsilon}^{\phi}$ as follows

[TABLE]

where $C$ is independent of $\epsilon$. Next, we use the result ([36, Lemma 5.1, p.7]), that is, $\left\|\phi_{0}\right\|_{H^{1}(U^{\epsilon})}\leq\epsilon^{1/2}C\left\|\nabla\phi_{0}\right\|_{H^{1}(\Omega^{\epsilon})}\,,$ for a $C$ independent of $\epsilon$. Herewith, we established (54). Using the trace theorem and the fact that \eta^{\phi}_{\epsilon}={\color[rgb]{0,0,0}{\rm G}^{\phi}_{\epsilon}} on $\partial\Omega^{\epsilon}$ allow us to obtain the estimate

[TABLE]

which provides with (54) the bound

[TABLE]

Applying the same arguments to ${\rm G}_{\epsilon}^{w}$ immediately leads to the corresponding bound

[TABLE]

Next, we estimate (41)1. Testing (41)1 with $-\Delta{\rm E}_{\epsilon}^{w}={\mathcal{A}}_{\epsilon}{\rm E}_{\epsilon}^{w}$ provides

[TABLE]

where we used (52) and the boundary terms $[BT1]$ and $[BT2]$ are bounded as follows

[TABLE]

In (60)1, we used Assumption C to assure in terms of regularity that $\phi_{0}=(-\Delta)^{-1}w_{0}\in C^{1}(0,T;W^{l,\infty}(\Omega^{\epsilon}))$ and hence the final bound is a consequence of the trace theorem and (58). The same argument holds for (60)2.

All this together then leads to the estimate

[TABLE]

where the last summand reflects the boundary terms. In order to control the terms on the right-hand side in (61), we make use of the fact that $f(s)$ is a polynomial and satisfies (50). The first term in (61) satisfies the following inequality

[TABLE]

Before we proceed, we estimate the terms on the right-hand side in (62):

1st term in (62): We first note that with the remainder term in Taylor series we obtain

[TABLE]

where $I_{\phi}:=]\underline{\phi},\overline{\phi}[$ is the interval defined by the smallest real root $\underline{\phi}$ and $\overline{\phi}$ the largest real root of the polynomial free energy $F$ characterized by (2). We used Assumption C in (63). Herewith, we can estimate the first term (e.g. in $d=3$) as follows

[TABLE]

2nd term in (62): With Sobolev inequalities, e.g. Fridrichs’ inequality in the perforated domain case [10], and the identity $\phi^{\epsilon}-\phi_{0}={\rm E}_{\epsilon}^{\phi}+\epsilon\phi_{1}$ we obtain the following estimate

[TABLE]

where we again used Assumption C and classical regularity results from the elliptic PDE theory.

3rd term in (62): Following the same ideas as for the 1st term estimated in (64), we immediately get the bound

[TABLE]

4th term in (62): The last term can finally be controlled as follows

[TABLE]

where we again used Assumption C and classical regularity results from the theory of elliptic PDEs.

Back to controlling (61), it leaves to control the second term on the right-hand side, i.e., $\left|\left({\cal A}_{\epsilon}[f^{\prime}(\phi_{0})({\rm E}^{\phi}_{\epsilon}+\epsilon\phi_{1})],{\cal A}_{\epsilon}E_{\epsilon}^{w}\right)\right|$. We have

[TABLE]

where we used the facts that

[TABLE]

Assumption C, and error equation (41)3.

Hence, with the previously derived bounds (62) and (68) we obtain,

[TABLE]

where $A(\epsilon):=\left(\epsilon^{2}(2+\epsilon^{1/2})^{2}+\epsilon^{2}+\epsilon+\epsilon^{1/2}\right)\,.$ A consideration of $\frac{d}{dt}\left({\rm exp}\left(-Ct\right)\left\|{\rm E}_{\epsilon}^{w}\right\|\right)$ leads after some rewriting to the following bound,

[TABLE]

which quarantees the control of (61). Herewith, we are also in the position to derive a bound on $\left\|{\cal A}_{\epsilon}^{w}E_{\epsilon}^{w}\right\|^{2}$ after integrating (70) over time, that means,

[TABLE]

for $n\in\mathbb{N}$ finite. $\square$

6 Conclusions

Based on a microscopic porous media formulation (8), we derived upscaled/ homogenized phase-field equations for general free energy densities. We gave a rigorous justification of this new effective macroscopic equations for a class of polynomial free energies which include the widely used double-well potential. The porous materials considered here can be represented by a periodic covering of a single reference cell $Y$ which accounts for the pore geometry. It is well-known that transport as well as fluid flow in porous media lead to high-dimensional computational problems, since the mesh size needs to be much smaller than the heterogeneity $\epsilon:=\frac{\ell}{L}$, the ratio of the characteristic length scale of the pores $\ell$ over the size of macroscopic porous medium $L$. We rigorously derived qualitative error estimates for the approximation error between the solution of the effective macroscopic problem (17) and the solution of the fully resolved microscopic equation (8). We also recovered the classical error behavior from homogenization of elliptic problems based on a truncation not resolving boundary effects in the context of fourth order phase field problems.

This error quantification is of fundamental interest in applications as it provides guidance on the applicability of the new effective macroscopic phase-field formulation in dependence of the heterogeneity $\epsilon>0$ defined by the heterogeneous material under consideration. Our dimensionally reduced phase field formulation can also be seen as the precursor to an effective and systematic computational strategy where microscopic properties such as the geometry and wall characteristics (e.g. wetting properties) of a reference pore enter the macroscopic description in an effective manner, thus avoiding a full numerical resolution of the finer details of the porous structure.

There are a number of interesting mathematical-physical questions related to the analysis presented here. For instance: (i) The error estimates still require rather high regularity assumptions (Assumption C) for which Lemma 1 seems to provide currently the best available estimates; (ii) Moreover, the error behavior in time in (20) seems not optimal; (iii) And a more physically motivated question is: “What is the influence of the pore or material geometry on phase transformations in heterogeneous media such as composites and porous materials and does the effective formulation capture such geometry dependent phenomena?”

Finally, we believe that, due to the popularity and the wide range of applicability of phase-field equations, the new effective macroscopic formulation, could ultimately serve as a promising computational tool in material-chemical- physical sciences and engineering. In particular, the effective phase-field formulation (17) could form the basis for a promising new direction for modelling multiphase flow in porous media without making use of Darcy’s law.

We shall examine these and related questions in future studies.

Acknowledgments

We would like to thank the Referees for their helpful comments and careful reading of the manuscript.