The Unique ergodic stationary distribution of two stochastic SEIVS epidemic models with higher order perturbation

1.   Introduction

The Cahn–Hilliard equation provides a continuous description of the phase separation process for binary mixtures. It was first proposed in ^[1,2,3] to model the so-called spinodal decomposition of binary alloys in a rapid cooling process, assuming isotropy of the material. As pointed out in ^[4], the Cahn–Hillirad equation is rather broad ranged in its evolution scope such that it is capable of describing important qualitative features of many systems undergoing phase separation at different time stages. Besides the spontaneous heterogenization of a binary mixture like spinodal decomposition, it can also model other mechanisms in pattern formation such as the process of nucleation and growth, and the process of coarsening ^[3,4,5]. The Cahn–Hilliard equation is a representative of the so-called diffuse interface models describing the evolution of free interfaces during phase transitions. Instead of the classical sharp-interface formulation, the diffuse interface model represents free interfaces that separate different components of the mixture as thin layers with finite thickness over which material properties vary smoothly. This approach has several advantages, see e.g., ^[6,7]. First, explicit tracking of the free interface, which is usually a difficult task, can be avoided in both mathematical formulation and numerical computation. Second, evolution of complex geometries and topological changes of the free interface can be handled in a natural way. In the past years, the Cahn–Hilliard equation and its variants have been successfully applied in different scientific fields and turned out to be efficient tools for the study of a wide variety of segregation-like phenomena ^[8], for instance, diblock copolymer ^[9], image inpainting ^[10,11], tumor growth ^[12,13], biology ^[14], two-phase flows ^{[7,15,16,17,18]}, moving contact line dynamics ^[19,20], and so on.

Assume that $T\in (0, +\infty]$, $\Omega\subset \mathbb{R}^d$ ($d = 1, 2, 3$) is a bounded domain with a sufficiently smooth boundary $\partial \Omega$. When $d = 1$, $\Omega$ is simply an interval, e.g., $\Omega = (0, l)$ for some $l > 0$. The classical Cahn–Hilliard equation can be written in the following form:

The unknown function $\phi$ is called the order parameter or phase-field, which is related to local concentrations of the two components of a binary mixture. It may have different interpretations according to the physical context, for instance, the volume fraction, mass fraction, or mole fraction ^[4,6]. Usually people consider a rescaled form of $\phi$ such that it denotes the difference between local concentrations for the two components, that is, $\phi = c_A-c_B$. Since the concentrations satisfy $c_A, c_B\in[0, 1]$ and $c_A+c_B = 1$, it is straightforward to check that $\phi$ should take its value in the physical interval $[-1, 1]$, with $\pm 1$ corresponding to the pure states. In a general framework for the description of free interfaces, the phase-field function $\phi$ takes distinct values in different bulk phases away from the diffuse interface separating them, and the free interface can be identified with an intermediate level set of $\phi$ (e.g., the zero-level set).

In Eq (1.1), $\mu$ denotes the chemical potential, $\epsilon$ is a positive constant and $M$ is a nonnegative quantity standing for the diffusion mobility that can be chosen as either a positive constant or a concentration dependent function $M = M(\phi)$. From the mathematical point of view, (1.1) is a fourth order parabolic equation for the unknown variable $\phi$ when $M > 0$. Thus, suitable initial and boundary conditions should be taken into account to form a well-posed problem. For the initial condition, we take

On the other hand, one of the classical choices for boundary conditions is the following homogeneous Neumann type:

Here, ${\bf{n}} = {\bf{n}}(x)$ denotes the unit exterior normal to the boundary $\partial\Omega$ and $\partial_{\bf{n}}$ stands for the outward normal derivative on $\partial \Omega$. Other boundary conditions that are of interest are Dirichlet boundary conditions (for $\phi$ and $\mu$) and periodic boundary conditions (cf. ^[21]).

The Cahn–Hilliard equation has been extensively studied in the literature both analytically and numerically. Our aim in this paper is to review the derivation, mathematical structure and some analytical issues for the Cahn–Hilliard equation and its variants, with main focus on the well-posedness and long-time behavior of global solutions. Some other important and interesting problems like nonlocal interactions, limit motions, formal and rigorous justifications of the sharp-interface limit etc are beyond the scope of this paper. In Section 2, we present results for the Cahn–Hilliard equation in the classical setting, i.e., the initial boundary value problem (1.1)–(1.4). In Section 3, we discuss some recent developments on the dynamic boundary conditions.

2.   Cahn–Hilliard equation in the classical setting

2.1.   Derivation and structure

For an isotropic binary mixture with nonuniform composition at a fixed temperature, we consider the following Ginzburg–Landau type free energy:

The energy functional $E^{\mathrm{bulk}}$ consists of two contributions, that is, the bulk part and the gradient part. The bulk energy part represents the interaction of different components in a homogeneous system, where $F(\phi)$ denotes the (Helmholtz) free energy density of mixing. A typical thermodynamically relevant example is the following logarithmic potential (see ^[2,4])

for some constants $\theta, \theta_c > 0$, where $\theta$ is the absolute temperature of the system and $\theta_c$ is the critical temperature of phase separation. The potential (2.2) is also related to the Flory–Huggins free energy density in assessing the mutual miscibility of polymer solutions. When $0 < \theta < \theta_c$, it is easy to verify that $F_{\mathrm{log}}$ has a double-well structure with two minima $\pm \phi_*\in (-1, 1)$, where $\phi_*$ is the positive root of the equation $F_{\mathrm{log}}'(s) = 0$ (see ^[22]). This case is of physical importance since the phase separation may occur, see the phase diagram in ^[4,23]. The interval $(-\phi^*, \phi^*)$ with $\phi^* = (1-\theta/\theta_c)^{1/2} > 0$ is called the spinodal region, in which it holds $F_{\mathrm{log}}''(s) < 0$. When a homogeneous state is located in $(-\phi^*, \phi^*)$, any small fluctuation in composition will lower the free energy and yield spontaneous phase separation of the mixture towards the equilibrium compositions of two phases corresponding to $\pm \phi_*$. When the system temperature $\theta$ is closed to the critical temperature $\theta_c$, i.e., the so-called "shallow quenching" case, the singular potential $F_{\mathrm{log}}$ is often approximated by a polynomial of degree four like the following form

More precisely, in view of (2.2), we can apply Taylor's expansion at $s = 0$ to get (cf. ^[24])

For the special case $\theta = 3$, $\theta_c = 4$, we recover (2.3) after adding the resultant with a constant $1/4$. The above simple approximation brings great convenience in the mathematical analysis and numeric simulation for the Cahn–Hilliard equation. Nevertheless, we note that it leads to a shift of the location of minima, i.e., from $\pm \phi_*$ to $\pm1$. In the literature, there is another type of singular potential called the double obstacle potential, which can be obtained in the "deep-quench limit" of the logarithmic potential (2.2) by letting $\theta\to 0^+$ (see ^[22]). In this case, the spinodal region eventually expands to $(-1, 1)$ and the corresponding free energy density usually takes the following form (after a shift by adding the constant $\theta_c/2$):

The gradient part of $E^{\mathrm{bulk}}$ represents the spatial variation in composition of the mixture, where $\epsilon^2 > 0$ is a parameter often called the coefficient of gradient energy ^[2,21,23]. The small positive constant $\epsilon$ measures the capillary width of the mixture (i.e., thickness of the transition layer), see ^[6,7]. Formally speaking, this gradient term can be obtained by making expansion of a free energy density in a region of nonuniform composition

see e.g., ^[2,23] for details. As pointed out in ^[25], the gradient part accounts for heterogeneity of the mixture and serves as a penalty for the phase-field function $\phi$ having sudden changes with respect to the spatial variable $x$. It corresponds to the tendency of the mixture to prefer to be uniform in space and gives an approximation of the interfacial surface energy ^[6]. Here, this gradient term can also be regarded as an elliptic regularization against the possible backwards diffusion (noting that $F$ can be nonconvex in the spinodal region), in order to guarantee the system of partial differential equations to be well-posed (cf. ^[4,25]). We note that the idea of regularization by gradient terms dates back to the earlier work ^[26].

2.1.1.   Derivation from the mass balance law

There are several ways to derive (and to understand) the Cahn–Hilliard equation (1.1). We first present the derivation by using basic thermodynamics (cf. ^[2,23]). Recalling the definition of the order parameter $\phi$, we see that local concentrations of the two components of a binary mixture can be determined once the scalar function $\phi$ is known. Thus, dynamics of the composition can be predicted by a single evolution equation for $\phi$.

The conserved dynamics of a phase separation process is due to the generalized (non-Fickian) diffusion driven by gradients in the chemical potential $\mu$ (see ^[2], and also ^[21,23,27]). Phenomenologically, we consider a continuity equation

where the vector ${\bf{J}}$ denotes the mass flux. Equation (2.5) yields a differential (local) form for the law of mass balance. The natural and possibly the simplest boundary conditions are as follows (cf. (1.3)–(1.4))

The first boundary condition (2.6) is usually called the no-flux boundary condition such that after integration by parts, we can easily deduce from (2.5) the mass conservation of the system (at least in a formal manner by assuming that the solution is smooth):

Therefore, it holds

On the other hand, the homogeneous Neumann boundary condition (2.7) for $\phi$ also has its physical interpretation: the free interface between the two components intersects the solid wall (i.e., the boundary $\partial \Omega$) at a perfect static contact angle of $\pi/2$ (cf. ^[19,20]). Under this choice, the chemical potential $\mu$, that is defined as the variational derivative of the bulk free energy $E^{\mathrm{bulk}}$ with respect to $\phi$, is given by

where $F'$ denotes the derivative of the bulk potential $F$ with respect to $\phi$ (cf. ^[23,25]). The chemical potential may be viewed as a forcing term proportional to the local distance from equilibrium ^[23]. Then the right expression of the mass flux ${\bf{J}}$ should be chosen to fulfill the basic thermodynamics, that is, the evolution of $\phi$ must occur in such a way that the free energy $E^{\mathrm{bulk}}$ does not increase in time. For instance, we may postulate the following constitutive equation

Thus, combining (2.5), (2.9) and (2.10), we arrive at the Cahn–Hilliard equation (1.1). Besides, it is obvious that (2.6) together with (2.10) yields the boundary condition (1.3). An important consequence of (2.6), (2.7) and (2.10) is that the free energy $E^{\mathrm{bulk}}$ is indeed non-increasing in time, which can be seen from a direct calculation (cf. ^[21,28]):

In this sense, (2.7) is sometimes referred to as the variational boundary condition such that it guarantees the validity of the energy dissipation property (2.11) and thus the resulting system fulfills the requirement from the laws of thermodynamics.

When the Cahn–Hilliard equation (1.1) is subject to other types of boundary conditions for $\phi$ and $\mu$, for instance, the periodic boundary conditions or the (nonhomogeneous) Dirichlet boundary conditions, a similar energy dissipation law can still be derived, see e.g., ^[21].

2.1.2.   Gradient flow structure

Another possible derivation of (1.1) we would like to mention is based on the gradient flow approach, see e.g., ^[25], where the author treated the case with $M$ being a positive constant. Within this framework, the Cahn–Hilliard equation (1.1) can be derived by considering the constrained gradient dynamics for the free energy $E^{\mathrm{bulk}}$, subject to the mass conservation property (2.8). Namely, we seek a law of evolution in the form

for some constant mobility $M > 0$. To specify the meaning of the constraint gradient denoted by ${\mathrm{grad}}_0$, we introduce some function spaces (cf. ^[25,28]). For every $f\in H^1(\Omega)'$, we denote by $\overline{f}$ its generalized mean value over $\Omega$ such that $\overline{f} = |\Omega|^{-1}\langle f, 1\rangle_{(H^1)', H^1}$; if $f\in L^1(\Omega)$, then its mean is simply given by $\overline{f} = |\Omega|^{-1}\int_\Omega f(x) \, {\mathrm{d}} x$. Consider the realization of the minus Laplacian with homogeneous Neumann boundary condition $\mathcal{A}_N\in \mathcal{L}(H^1(\Omega), H^1(\Omega)')$ defined by

For the linear spaces

we can check that the operator $\mathcal{A}_N$ is self-adjoint, positively defined on $\dot{H}^1(\Omega)$ and the restriction of $\mathcal{A}_N$ from $\dot{H}^1(\Omega)$ onto $\dot{H}^{-1}(\Omega)$ is an isomorphism (an easy consequence from the Poincaré–Wirtinger inequality and the Lax–Milgram theorem). Denote its inverse by $\mathcal{N} = \mathcal{A}_N^{-1}: \dot{H}^{-1}(\Omega) \to \dot{H}^{1}(\Omega)$. Then for every $f\in \dot{H}^{-1}(\Omega)$, $u = \mathcal{N} f \in \dot{H}^{1}(\Omega)$ is the unique weak solution of the Neumann problem

Moreover, it holds

Let

We remark that the second boundary condition for $\Delta \phi$ is actually equivalent to $\partial_{\bf{n}}\mu = 0$ on $\partial\Omega$ thanks to (1.3)–(1.4). Let $\phi\in H^4_N(\Omega)$. Then for any smooth function $f$ satisfying $\overline{f} = 0$ and $\partial_{\bf{n}} f = 0$ on $\partial\Omega$, we use the expression $f = -\Delta(\mathcal{N} f)$ and integration by parts to obtain (in the sense of Gateaux derivative, cf. [28,Section 3])

Hence, we can identify in $\dot{H}^{-1}(\Omega)$ that

and specify its domain as $H^4_N(\Omega)$. From (2.12) and (2.13) we obtain

which is exactly the Cahn–Hilliard equation (1.1) with the mobility $M$ being a positive constant and subject to the boundary conditions (1.3)–(1.4).

The above simple argument indicates that in case of a constant mobility, the Cahn–Hilliard equation can be regarded as a constrained gradient flow in the Hilbert space $\dot{H}^{-1}(\Omega)$. Concerning the more general case with concentration dependent mobilities, it was shown in ^[29] that for the linear mobility $M(s) = s$, the Cahn–Hilliard equation generates a gradient flow in the space endowed with a non-Hilbertian metric, i.e., the $L^2$-Wasserstein distance (note that the solution considered in ^[29] is nonnegative in the framework therein). Moreover, in ^[30], the authors considered the Cahn–Hilliard type equation with some general nonlinear mobilities (e.g., nonnegative concave functions) as a gradient flow in certain weighted-Wasserstein metric spaces and proved the existence of weak solutions by the variational minimizing movement approach.

2.1.3.   Derivation by the energetic variational approach

We note that the previous derivations, though performed along different procedures, aim to derive a set of partial differential equations that mainly fulfill two physical constraints: the mass conservation (expressed either in a local or a global form) and the energy dissipation. The constitutive relation (2.10) is postulated to guarantee the energy dissipation (2.11), which may not be the unique choice for this purpose. There are some other evolution equations that have similar properties, for instance, the conservative Allen–Cahn equation (see ^[31])

subject to the homogeneous Neumann boundary condition $\partial_{\bf{n}}\phi = 0$ on $\partial \Omega\times (0, T)$. In this case, a direct calculation yields that (2.8) holds and the free energy $E^{\mathrm{bulk}}$ is still non-increasing in time:

Below we provide an alternative derivation of the Cahn–Hilliard equation (1.1) via the energetic variational approach, which combines the least action principle and Onsager's principle of maximum energy dissipation in continuum mechanics ^[32,33,34]. Within this general framework, one can easily include different physical processes by choosing specific form of the free energy as well as the dissipation for the nonequilibrium system. Then based on suitable kinematic and energetic assumptions, the corresponding partial differential equations can be uniquely determined by force balance relations. This approach has been successfully applied to derive complex hydrodynamic systems in fluid dynamics, liquid crystals, electrokinetics, visco-elasticity, and so on, we refer to ^{[35,36,37,38,39]} and the references cited therein.

In the domain $\Omega$, $\phi$ is assumed to be a locally conserved quantity that satisfies the continuity equation

Here, the mass flux is given by ${\bf{J}} = \phi {\bf{u}}$, where ${\bf{u}}: \Omega \to \mathbb{R}^d$ stands for the microscopic effective velocity (e.g., due to certain diffusion process). Concerning the boundary conditions, we again assume the no-flux boundary condition ${\bf{J}}\cdot {\bf{n}} = 0$ on $\partial \Omega$, which can be guaranteed by

For an isothermal closed system, we assume that the evolution of a binary mixture satisfies the following energy dissipation law, which is a natural consequence of the first and second laws of thermodynamics (see e.g., ^[36]):

where the bulk free energy $E^{{\mathrm{bulk}}}$ takes the form in (2.1) and the rate of energy dissipation $\mathcal{D}^{{\mathrm{bulk}}}$ is chosen as

with $M > 0$ being the mobility.

It remains to determine the microscopic velocity ${\bf{u}}$ in (2.15). Let $\Omega_0^X, \Omega_t^x\subset\mathbb{R}^d$ be bounded domains with smooth boundaries $\Gamma_0^X, \Gamma_t^x$, respectively. We introduce the flow map $x(X, t): \Omega_0^X\to \Omega_t^x$, which is defined as a solution to the system of ordinary differential equations

where $X = (X_1, ..., X_d)^\mathrm{tr}\in \Omega_0^X$, $x = (x_1, ..., x_d)^\mathrm{tr}\in \Omega_t^x$, and ${\bf{u}}(x, t)\in \mathbb{R}^d$ is a (sufficiently smooth) velocity field. The coordinate system $X$ is called the Lagrangian coordinate system and it refers to $\Omega_0^X$ that we call the reference configuration. Meanwhile, the coordinate system $x$ is called the Eulerian coordinate system and it refers to $\Omega_x^t$ that we call the deformed configuration. We shall denote $\nabla_x$ the spatial gradient operator in $\Omega$ under the Eulerian coordinate system. Introduce now the action functional

Then applying the least action principle, which states that the dynamics of a Hamiltonian system is determined by a critical point of the action functional with respect to the trajectory (in the Lagrangian coordinates), we eventually get

In the above computation, we have assumed a trivial boundary dynamics for $\phi$ such that the boundary condition (1.4) holds. The relation (2.21) yields the conservative force (written in the strong form, cf. [38,Section 2.2.1])

On the other hand, from Onsager's maximum dissipation principle for a dissipative system, we can calculate the dissipative force by taking variation of the Rayleigh dissipation functional $\mathcal{R} = \frac12\mathcal{D}^{{\mathrm{bulk}}}$ (recall (2.18)) with respect to the rate function ${\bf{u}}$, that is,

This gives the generalized dissipative force (again written in the strong form, cf. [38,Section 2.2.1]):

By the force balance relation, i.e., Newton's second law ${\bf{f}}_{\rm{inertial}}+{\bf{f}}_{\rm{conv}}+{\bf{f}}_{\rm{diss}} = 0$ (recalling that here we have ${\bf{f}}_{\rm{inertial}} = 0$ because the kinetic energy is neglected), we obtain

where $\mu = -\epsilon^2 \Delta \phi +F'(\phi)$. Finally, solving ${\bf{u}}$ from the above algebraic equation and inserting it back into (2.15), we arrive at the Cahn–Hilliard equation (1.1), subject to the boundary conditions (1.3)–(1.4).

The above derivation via the energetic variational approach reveals that the Cahn–Hilliard equation together with the classical boundary conditions naturally fulfills three important physical constraints: conservation of mass, dissipation of energy and, in addition, the force balance.

2.2.   Well-posedness

The Cahn–Hilliard equation (1.1) subject to the initial and boundary conditions (1.2)–(1.4) has been extensively studied in the literature. A rather complete picture about the existence, uniqueness, regularity and long-time behavior of global solutions has been obtained. For details, we refer the reader to the recent book ^[41] and the references cited therein.

2.2.1.   The case with a regular potential

Concerning problem (1.1)–(1.4) with a constant mobility (e.g., $M = 1$ without loss of generality) and a regular potential $F$ (e.g., a fourth order polynomial like (2.3)), existence of global weak solutions with an initial datum $\phi_0\in H^1(\Omega)$ can be easily obtained by using the Faedo–Galerkin method, thanks to the Lyapunov structure (2.11) (see ^[21,42]). Existence and uniqueness of more regular solutions with an initial datum $\phi_0\in H^2_N(\Omega) = \big\{\phi\in H^2(\Omega)\, :\ \partial_{\bf{n}}\phi = 0\ {\rm{on}}\ \Omega\big\}$ was proved in ^[91]. More precisely, they considered

where $\gamma_1, \gamma_2$ are two constant parameters. It was shown that the sign of $\gamma_2$ (coefficient of the leading order term) plays a crucial role in the study of existence of global solutions. If $\gamma_2 > 0$, for any initial datum there is a unique global solution, while for $\gamma_2 < 0$, the solution must blow up in a finite time for large initial data. The proofs therein relied on the Picard iteration scheme (for local existence and uniqueness) and the energy method (for global existence and finite time blow-up). Besides, the authors studied a finite element Galerkin approximation of the initial boundary value problem and obtained the existence result as well as some optimal order error bounds. We also refer to ^[44], in which the authors considered a general polynomial of the following form

and boundary conditions of either Neumann or periodic type. The specific form and growth condition on the regular potential $F$ can be relaxed when studying the regularity of solutions, see for instance, [41,Chapter 3,Section 3.4] such that $F\in C^3(\mathbb{R})$, $F(0) = F'(0) = 0$ and

Besides, in ^[45], the author analyzed the Cahn–Hilliard equation (1.1) with $M = 1$, subject to (1.2) and the homogeneous Dirichlet boundary conditions $\phi = \Delta\phi = 0$ on $\partial\Omega$, where $\Omega \subset\mathbb{R}^d$, $d = 1, 2, 3$, is a bounded domain with smooth boundary. He proved the unique global solvability in $H^1_0(\Omega)$ and the existence of a global attractor, for a general regular function $F\in C^3(\mathbb{R})$ satisfying $F(0) = 0$ and

where $C\geq 0$, $p\leq 1$ when $d = 3$ and $p$ is arbitrary when $d = 2$. In the above assumptions, the growth condition on the nonlinearity is essentially due to the Sobolev embedding theorem.

It is worth noting that when the mobility $M$ is a positive constant, (1.1) is a fourth-order semilinear parabolic equation. In general, a fourth-order parabolic equation does not maintain the maximum principle that holds for second-order parabolic equations. Concerning the equation (1.1), its solution $\phi$ needs not to stay in the physically relevant interval $[-1, 1]$ as time evolves, even if the initial datum $\phi_0$ satisfies this property. We can refer to [41,Chapter 4,Remark 4.10] for a simple counterexample regarding the Cahn–Hilliard equation with $M = 1$ and $F'(\phi) = \phi^3-\phi$ in the case of one space dimension.

Next, we would like to mention the viscous Cahn–Hilliard equation

with the viscous parameter $\alpha > 0$. Equation (2.22) was introduced in ^[46] to include certain viscoelastic relaxation effects in the phase separation process, which was neglected in the original work ^[2]. The viscous term $\alpha \partial_t \phi$ is also related to the influence of certain internal microforces in the mixture (see ^[47]). From the energetic point of view, it yields some additional dissipation in the system, for instance, any smooth solution to the initial boundary value problem (2.22) with (1.2)–(1.4) satisfies

Thus, the viscous Cahn–Hilliard equation (2.22) can be regarded as a regularized version of the original Cahn–Hilliard equation (1.1). When the potential $F$ is regular, e.g., a general polynomial, we refer to ^[43,48] for an extensive study on the computation and mathematical analysis for the equation (2.22) subject to the homogeneous Dirichlet (or Neumann) boundary conditions as well as the initial condition (1.2). The results in ^[43,48] showed that in a suitable sense the viscous Cahn–Hilliard equation (2.22) can actually be viewed as an interpolation between the Cahn–Hilliard equation (1.1) (with $M = 1$) and the Allen–Cahn equation for the grain-boundary migration (see ^[40]), that is,

2.2.2.   The case with a degenerate mobility

Next, we consider the situation that the diffusion mobility $M$ is concentration dependent, and even may be degenerate at some values of $\phi$. This case appeared in the original derivation of the Cahn–Hilliard equation ^[1] and a thermodynamically reasonable choice is $M(s) = 1-s^2$ (see ^[4,23,67,68]). Concerning the mathematical analysis, the author of ^[49] studied a general Cahn–Hilliard type equation in one space dimension (i.e., on the interval $(0, 1)$) with $F'$ being continuous and $M$ being Hölder continuous such that

Then for any initial datum $\phi_0\in H^3_0(0, 1) = \big\{\phi\in H^3(0, 1)\, :\, Du = 0\ \rm{at}\ x = 0, 1\big\}$ satisfying $0\leq \phi_0\leq 1$, he proved the existence of a global generalized solution $\phi$ with a uniform bound in $L^\infty(0, T;H^1(0, 1))$ and a local $L^2$-estimate on $D^3\phi$ when the equation does not degenerate. In particular, thanks to the degeneracy of the mobility, he obtained nonnegativity of the solution such that $0\leq \phi(x, t)\leq 1$, see [49,Theorem 3.2]. We note that in ^[49], the phase function $\phi$ stands for the concentration of one component, not the difference of concentrations. The result of ^[49] implies that the degenerate mobility $M$ turns out to be a sufficient condition for the existence of weak solutions with the physical property $0\leq \phi(x, t)\leq 1$ for $(x, t)\in (0, 1)\times (0, +\infty)$, as long as $0\leq \phi_0(x)\leq 1$ for $x\in (0, 1)$.

Later on, the Cahn–Hilliard equation (1.1) with a degenerate concentration dependent mobility as well as a singular potential was analyzed in ^[50]. Assume that

(ⅰ) $F(s) = F_1(s)+F_2(s)$, with $\|F_2(s)\|_{C^2[-1, 1]}\leq C$ and $F_1:(-1, 1)\to \mathbb{R}$ being convex, satisfying $F''_1(s) = (1-u^2)^{-m}\widetilde{F}(s)$ for $m\geq 1$ and a $C^1$-function $\widetilde{F}:[-1, 1]\to \mathbb{R}^+\cup\{0\}$;

(ⅱ) $M(s) = (1-s^2)^m\widetilde{M}(s)$ with a $C^1$-function $\widetilde{M}:[-1, 1]\to \mathbb{R}$ satisfying $0 < m_0\leq \widetilde{M}(s)\leq M_0$ for $s\in [-1, 1]$, and $M(s)$ is extended to $\mathbb{R}$ by setting $M(s) = 0$ for $|s| > 1$.

It is easy to check that these general assumptions cover the physically relevant case with the logarithmic potential (2.2) and $M(s) = 1-s^2$, $s\in[-1, 1]$. Let either $\partial \Omega\in C^{1, 1}$ or $\Omega$ be convex in $\mathbb{R}^d$ ($d\geq 1$). Then for any initial datum $\phi_0\in H^1(\Omega)$ with $|\phi_0|\leq 1$ almost everywhere in $\Omega$ and $\int_\Omega F(\phi_0)+ \Phi(\phi_0)\, \mathrm{d}x < +\infty$, problem (1.1)–(1.4) admits a global weak solution on an arbitrary time interval $[0, T]$, see [50,Theorem 1]. Here, the function $\Phi:(-1, 1)\to \mathbb{R}^+\cup \{0\}$ is determined by

The proof therein relied on some suitable approximation of the degenerate mobility and a specific regularization of the singular potential such that (see [50,Section 3])

$F_\varepsilon = F^1_\varepsilon+F_2$ with

and $F_2$ being extended to be a function on $\mathbb{R}$ such that $\|F_2\|_{C^2(\mathbb{R})}\leq C$. We note that the special structure that $M(s)F''(s)$ is bounded plays a crucial role in the proof of existence in ^[50].

Recently, the authors of ^[51] considered the Gibbs–Thomson effect in the phase separation process. For the Cahn–Hilliard equation (1.1) in $\Omega = [0, 2\pi]^d$ subject to periodic boundary conditions, assuming that the double-well potential $F$ is smooth at its minima $\pm1$ (cf. (2.3)) and the mobility is of the form $M(s) = |1-s^2|^m$ without zero extension outside $[-1, 1]$ ($m > 0$ when $d = 1, 2$, and $m\in (0, 2/(d-2))$ when $d\geq 3$), they proved that for any initial datum $\phi_0\in H^1(\Omega)$, the resulting initial boundary value problem admits a global weak solution on an arbitrary time interval $[0, T]$, see [51,Theorem 2]. The proof is again based on some regularization of the degenerate mobility and then passing to the limit in the non-degenerate approximating problem. One difference from the previous results (cf. those in ^[50]) is that for the solution $\phi$, even if its initial value lies in $[-1, 1]$, it may not remain inside $[-1, 1]$ as long as the interface has nonzero mean curvature, which accommodates the physical Gibbs–Thomson effect. Nevertheless, we remark that in all the studies mentioned above, the uniqueness of weak solutions (in various formulations) to the Cahn–Hilliard equation with a degenerate mobility remains an open problem.

2.2.3.   The case with a singular potential

A lot of attentions have been paid to the problem (1.1)–(1.4) with a constant mobility $M > 0$ and a singular potential $F$ including the typical ones (2.2) and (2.4), see ^[27,41,58] for some surveys on this topic.

To study the existence of solutions, one basic strategy is to regularize the singular potential in a suitable manner, prove the existence of solutions to the regularized problem, derive uniform estimates with respect to the approximating parameter and then pass to the limit to extract a convergent subsequence.

In ^[22], the authors studied problem (1.1)–(1.4) in a bounded smooth domain $\Omega\subset \mathbb{R}^d$ ($d = 1, 2, 3$) with the double obstacle potential (2.4). They proved that for any initial datum $\phi_0\in H^1(\Omega)$, $|\phi_0|\leq 1$ almost everywhere in $\Omega$ with its spatial average in $(-1, 1)$, problem (1.1)–(1.4) admits a unique global weak solution $\phi$ satisfying $-1\leq \phi(x, t)\leq 1$ almost everywhere in $\Omega \times (0, T)$, see [22,Theorem 2.2]. Some regularity results on the weak solutions were also obtained in [22,Section 2.2]. The proof therein relied on the following $C^2$-regularization for the singular potential (see [22,(2.9)]) such that for $\varepsilon \in (0, 1)$,

When the logarithmic potential (2.2) is concerned, several different methods have been developed in the literature. The authors of ^[52] studied a general multi-component problem in a bounded smooth domain $\Omega\subset \mathbb{R}^d$ ($d = 1, 2, 3$) subject to homogeneous Neumann boundary conditions. Their results can easily apply to the two-component problem, that is, for potentials like

if the initial datum $\phi_0\in H^1(\Omega)$, $0\leq \phi_0(x)\leq 1$ almost everywhere in $\Omega$ with its spatial average belonging to $(0, 1)$, problem (1.1)–(1.4) admits a unique global weak solution $\phi$ satisfying $0\leq \phi(x, t)\leq 1$ almost everywhere in $\Omega \times (0, T)$, see [52,Theorem 1] for a general statement of the result. Indeed, the arguments in [52,Section 3] also implied that $0 < \phi(x, t) < 1$ almost everywhere in $\Omega \times (0, T)$, i.e., the set of singular points for $F'(s)$ has zero measure. Besides, in [52,Theorem 2] they justified the "deep-quench limit" problem studied in ^[22] by letting $\theta\to 0^+$. The key idea of the proof in ^[52] is to regularize the function $G(s) = s\ln s$ by

for some $\varepsilon\in (0, 1)$. After solving the regularized problem, they derived estimates that are uniform in $\varepsilon$ and then passed to the limit as $\varepsilon\to 0^+$ to extract a convergent subsequence, whose limit is a global weak solution to the original problem.

Later in ^[24], the authors introduced a different approximation of the logarithmic potential (2.2) in terms of the following polynomials

Then for any initial datum $\phi_0\in L^2(\Omega)$ (or $H^1(\Omega)$), $\|\phi_0\|_{L^\infty(\Omega)}\leq 1$ with its spatial average in $(-1, 1)$, they prove that the initial boundary value problem (1.1)–(1.4) (or in a periodic setting) admits a unique global weak solution satisfying $\|\phi(t)\|_{L^\infty(\Omega)}\leq 1$ for $t\geq 0$ and for $t > 0$, the set $\{x\in \Omega\, :\ |\phi(x, t)| = 1\}$ has measure zero, see [24,Theorem 2.2]. Moreover, they showed the continuity property of weak solutions, which yields the existence of a $C_0$-semigroup defined by $S(t):\phi_0\to \phi(t)$ in suitable phase spaces.

In ^[27], the author considered the following equation

with some constant $\chi\geq 0$ and a general class of singular nonlinearities $F'\in C^1(-1, 1)$ satisfying $F(0) = F'(0) = 0$ and

In particular, the logarithmic potential (2.2) fulfills the above assumptions. By adapting the methods in ^[50,52], he made use of the following approximation $F'_\varepsilon\in C^1(\mathbb{R})$ such that for $\varepsilon\in (0, 1)$,

Let $\Omega\subset \mathbb{R}^d$, $d = 1, 2, 3$ be a smooth bounded domain. The author proved that for any initial datum $\phi_0\in H^1(\Omega)$, $|\phi_0(x)| < 1$ almost everywhere in $\Omega$ with its spatial average belonging to $(-1, 1)$, the initial boundary value problem (2.24) with (1.2)–(1.4) admits a unique global weak solution $\phi$ satisfying $-1 < \phi(x, t) < 1$ a.e. in $\Omega \times (0, T)$, see [27,Theorem 2.6,Remark 2.7,Remark 2.8]. We note that when $\chi > 0$, the equation (2.24) is usually referred to as the Cahn–Hilliard–Oono equation, which was introduced in ^[53] to model certain long-ranged (nonlocal) interactions. For further details on its mathematical analysis, see e.g., ^[54] (with the regular potential (2.3)) and ^[55] (with singular potentials including (2.2)). In ^[54], it was also shown that, in certain sense (i.e., the existence of a robust family of exponential attractors as $\chi\to0^+$), the dynamics of the original Cahn–Hilliard equation is "close" to that of the Cahn–Hilliard–Oono equation, for $\chi > 0$ being sufficiently small.

Besides the approximating arguments described above, there are some other approaches to handle the Cahn–Hilliard equation with singular potentials. We first mention the work ^[56], which treated problem (1.1)–(1.4) (again with a constant mobility) in a rather general setting. For a wide class of non-smooth potentials including both (2.2) and (2.4), thanks to the decomposition

where $\widehat{\beta}$ is a proper, lower semi-continuous and convex function and $\widehat{g}$ is a $C^1$ function with Lipschitz continuous derivative $g = \widehat{g}'$, the Cahn–Hilliard equation (1.1) can be written in a generalized form involving multi-valued mappings via subdifferential operators. In this framework, the authors of ^[56] first regularized the Cahn–Hilliard equation by the viscous one (cf. (2.22)) written in the following form

where $\alpha > 0$ and $\mathcal{N}$ denotes the inverse of minus Laplacian subject to the homogeneous Neumann boundary condition on the space $L^2_0(\Omega) = \big\{\phi\in L^2(\Omega): \int_\Omega \phi(x)\, {\mathrm{d}} x = 0\big\}$. Then they solved the viscous problem (indeed in a more general form) using an argument based on Yosida's approximation of monotone operators. After deriving uniform estimates for the approximate solutions with respect to the viscous parameter $\alpha$, they obtained the existence of a global weak solution to the original problem (1.1)–(1.4) by finding a convergent subsequence as $\alpha\to 0^+$ (see [56,Theorem 6.1] for a complete statement including also results on uniqueness and long-time behavior). On the other hand, with the help of monotone operator methods ^[57], the authors of ^[65] provided a different proof for the existence and uniqueness of global weak solutions to problem (1.1)–(1.4) with logarithmic type potentials (see [65,Theorem 1.2]). After decomposing the nonlinear singular term $F'$ in (1.1) into a monotone operator plus a globally Lipschitz continuous one, they achieved the conclusion in a direct manner by solving an abstract Cauchy problem for a suitable Lipschitz perturbation of a monotone operator, which is the subgradient of the convex part of the energy $E^{{\mathrm{bulk}}}$ (see [65,Theorem 3.1]).

2.2.4.   Regularity of solutions and separation from the pure states

Regularity of solutions to problem (1.1)–(1.4) has been investigated in some works mentioned above and the references cited therein. Roughly speaking, the parabolic nature of the Cahn–Hilliard equation (valid when $M > 0$) can yield some instantaneous regularizing effect of its weak solutions for $t > 0$. This fact is well understood for the case with regular potentials, while the case with singular potentials like (2.2) and (2.4) turns out to be more tricky. On one hand, the possible singularity of the potential (or its derivatives) at the pure states $\pm 1$ guarantees that the order parameter $\phi$ stays in the physically relevant interval $(-1, 1)$ (cf. the case with a degenerate mobility). On the other hand, it brings further difficulties to gain higher-order spatial regularity of the solution, see extensive discussions made in ^{[22,24,27,41,56,58,59,65]}.

The case with a logarithmic potential like (2.2) is of particular interest, because its two minima $\pm\phi_*$, which are the two nonzero solutions to the equation

locate exactly inside the interval $(-1, 1)$. From the dissipative nature of the Cahn–Hilliard equation, if the initial value is not a pure state (e.g., $|\Omega|^{-1}|\int_\Omega \phi_0(x)\, {\mathrm{d}} x| < 1$), one may expect in this case the following strict separation property along evolution:

for some constant $\delta\in (0, 1)$. The property (2.26), if it holds, implies that in the conserved phase separation process governed by the Cahn–Hilliard equation, the pure states can never be completely reached. From the mathematical point of view, it plays an important role in the analysis of (1.1), since the singular potential (2.2) can thus be regarded as a smooth, globally Lipschitz function so that further regularity of solutions to (1.1) can be obtained, see e.g., ^{[27,41,59,65]}.

In ^[59], the authors considered the viscous Cahn–Hilliard equation (2.22) when the spatial dimension is less or equal than three. Taking advantage of the viscous term $\alpha \partial_t \phi$ (with $\alpha > 0$) and using the comparison principle for second-order parabolic equations, they proved that the separation property (2.26) holds for all $t > 0$ (i.e., the so-called instantaneous strict separation property) and the separation is indeed uniform when $t\geq \eta$, for an arbitrary but fixed $\eta > 0$ (see [59,Corollary 3.2]). However, the distance $\delta$ from the pure states $\pm 1$ depends on the size of the viscous parameter $\alpha$ so that one cannot pass to the limit as $\alpha\to 0^+$ to draw the same conclusion for the original Cahn–Hilliard equation (1.1).

When $\alpha = 0$, the spatial dimension indeed plays a crucial role in the argument. In one and two dimensions, the instantaneous strict separation property was proved in [59,Proposition 7.1,Theorem 7.2] for general singular potentials satisfying the conditions in (2.25). When the spatial dimension is two, an additional assumption was also required therein:

with some positive constants $C_1$, $C_2$. The above condition allows one to apply the Trudinger–Moser type inequality (see e.g., ^[60]). Nevertheless, the assumptions mentioned above are fulfilled by the logarithmic potential (2.2). Later on, the authors of ^[55] extended the result to the two dimensional Cahn–Hilliard–Oono equation (2.24) by employing an alternative approach, which relies on the analysis of the following Neumann problem of an elliptic equation with a singular nonlinearity:

Here, the nonlinear function $\widetilde{F}$ is the strictly convex part of $F$ satisfying $\lim_{s\to \pm 1} \widetilde{F}(s) = \pm \infty$ (see [55,Section 3]). We refer to ^{[61,62,63,64]} for generalizations to the higher order Cahn–Hilliard type equation and some coupled systems with fluid interactions. In ^[55], an essential step to derive the separation property of $\phi$ is to prove that the singular derivative satisfies $\widetilde{F}'(\phi)\in L^\infty(\Omega\times(0, T))$, which follows from the fact $\mu\in L^\infty(0, T; H^2(\Omega))$ (this corresponds to assuming $f\in H^2(\Omega)$ in the elliptic problem (2.27)). In the recent work ^[63], the requirement on the $L^\infty_t H_x^2$-regularity of the chemical potential $\mu$ was weakened to be $\mu\in L^\infty(0, T; H^1(\Omega))$ by applying a suitable chain rule concerning the singular nonlinearity $\widetilde{F}'(\phi)$ (see [63,Lemma 3.2]).

When the spatial dimension is three, the situation is less satisfactory, because the singularity of the logarithmic potential at $\pm 1$ seems not strong enough. In ^[59], the authors obtained the instantaneous strict separation property, under a stronger assumption on the potential such that (2.25) is satisfied together with

for some $C > 0$ (see also ^[127] for some recent improvements). The above assumption is valid, for instance, for a class of nonlinearities like (see [59,Remark 7.1])

but it is not satisfied by the logarithmic potential (2.2) which is physically important. On the other hand, in ^[65], using the strict separation property for solutions to the stationary problem of the Cahn–Hilliard equation (see [65,Proposition 6.1]), the authors proved that global weak solutions of problem (1.1)–(1.4) with logarithmic type potential will eventually stay away from the pure states $\pm 1$ for sufficiently large time. Besides, the stability result for the Cahn–Hilliard–Hele–Shaw system that was obtained in [61,Theorem 2.5] implies that (after simply neglecting the fluid interaction), if the initial datum $\phi_0$ is not a pure state and belongs to a suitably small $H^2$-neighborhood of a local minimizer of the free energy $E^{{\mathrm{bulk}}}$, then the unique global strong solution to problem (1.1)–(1.4) with a logarithmic potential exists and is strictly separated from the pure phases $\pm 1$ for all $t\geq 0$.

2.3.   Long-time behavior

In what follows, we review some results on the long-time behavior of global solutions to problem (1.1)–(1.4) as $t\to +\infty$. Generally speaking, the study of long-time behavior of global solutions to a nonlinear evolution equation can be divided into two categories: the first category is to investigate the long-time behavior of the solution corresponding to a given initial datum, while the second category is to investigate the long-time behavior of a bundle of solutions whose initial data vary in a bounded set (see ^[146]). The latter category is related to the study of the associated infinite dimensional dynamical system ^[140], for instance, the existence of (finite dimensional) global attractors, exponential attractors, inertial manifolds, and so on. In this direction, there is a huge literature on problem (1.1)–(1.4) and its variants, we refer to ^{[24,27,42,44,45,54,55,58,59,101,136,138,140]} and the references cited therein. In particular, we refer to the recent book ^[41] for more detailed information.

In this paper, we confine ourselves to the discussion on the long-time behavior of a single trajectory defined by the global solution to problem (1.1)–(1.4). In view of the energy dissipation relation (2.11), it is natural to ask whether every bounded global solution to the evolution problem (1.1)–(1.4) will converge to a single equilibrium (e.g., local or global minimizers of the energy functional $E^{{\mathrm{bulk}}}$) as $t\to+\infty$. This property is sometimes called the uniqueness of asymptotic limit as time tends to infinity. In the one dimensional case, the long-time behavior of global solutions to problem (1.1)–(1.4) with a regular potential like (2.3) was first analyzed in ^[145]. Based on the well-posedness result obtained in ^[91] and the Lyapunov structure (2.11), the author proved that for any initial datum $\phi_0\in H^4_N(\Omega)$, the trajectory of solution is relatively compact in $H^2(\Omega)$ and moreover, when $t\to+\infty$, $\phi(x, t)$ converges to the $\omega$-limit set of $\phi_0$ denoted by $\omega(\phi_0)$, which is a compact, connected subset of $H^2(\Omega)$ and is positive invariant under the nonlinear semigroup $S(t):\phi_0 \to \phi(t)$ defined by the solution (see [145,Theorem 2.2]). Besides, he proved that each element of $\omega(\phi_0)$ is an equilibrium and in the one dimensional case, the associated stationary problem has only finite number of solutions (see [145,Theorem 3.4]). As a consequence, the connected set $\omega(\phi_0)$ consists of only one point so that the global solution $\phi(t)$ to the time-dependent Cahn–Hilliard equation must converge to an equilibrium as $t \to +\infty$ (see [145,Theorem 3.5]). In this regard, we also refer to ^[137] for an extended result for the one dimensional non-isothermal Cahn–Hilliard system.

The situation in higher spatial dimensions is more complicated. Although the dissipative structure of problem (1.1)–(1.4) (cf. (2.11)) guarantees that the $\omega$-limit set of a given initial datum $\phi_0$ contains only the steady states. However, this property is not sufficient to show that $\omega(\phi_0)$ is a singleton. One difficulty is that, in higher dimensions, the study on stationary solutions to the Cahn–Hilliard equation turns out to be rather involved (see e.g., ^[72,141]). In particular, due to the possible non-convexity of the free energy $E^{{\mathrm{bulk}}}$, the structure of the set of equilibria seems far from being well understood and estimates on the number of equilibria are only available in the one dimensional case so far ^[118,145]. On the other hand, examples have been given in the literature such that even for certain simple semilinear second order parabolic equation with a specific smooth nonlinear term, it admits a globally bounded solution whose $\omega$-limit set is diffeomorphic to the unit circle $S^1$, namely, a continuum (see e.g., ^[131]).

Various assumptions have been made in the literature to assure that a bounded global solution of an evolution equation will converge to a single equilibrium as $t\to +\infty$. Among these attempts, an efficient approach was proposed by ^[139], in which the author generalized the Łojasiewicz inequality for analytic functions in the finite-dimensional space $\mathbb{R}^m$ to the infinite dimensional spaces and proved that, if the nonlinear term is real analytic in the unknown function then the convergence to a single equilibrium as $t\to +\infty$ holds.

For problem (1.1)–(1.4) with a constant mobility ($M = 1$) and a general regular potential $F:\Omega\times \mathbb{R}\to \mathbb{R}$ depending analytically on $x$, $\phi$ uniformly in $x$ with suitable growth assumptions, the authors of ^[135] applied the Łojasiewicz–Simon approach and proved that the unique global solution converges to a single equilibrium as $t\to +\infty$ in the topology of $C^{4, \mu}(\Omega)$ with $\mu\in (0, 1)$ (see [135,Theorem 3.1]). The proof therein relied on a generalized gradient inequality of Łojasiewicz–Simon type (see [135,Theorem 3.2]), which is applicable to the $H^{-1}$-gradient flow defined by the Cahn–Hilliard equation (cf. (2.12)). Indeed, from the energy dissipation relation (2.11) (with $M = 1$) and the equation (1.1), we see that

On the other hand, from the precompactness of the trajectory $\phi(t)$ in certain space, for instance, $C^{4, \mu}(\Omega)$ as shown in ^[135], one can find a convergent subsequence $\phi(t_n)\to \phi_\infty$ in $C^{4, \mu}(\Omega)$, as $t_n\nearrow +\infty$ for some equilibrium $\phi_\infty\in \omega(\phi_0)$ (thanks to the characterization on $\omega(\phi_0)$ for gradient systems). At this stage, if one can show (which is not obvious indeed, since $h(t)\in L^2(0, +\infty)$ does not imply $h(t)\in L^1(0, +\infty)$)

then it holds

As a consequence, we have

which further implies the convergence of $\phi(t)$ to $\phi_\infty$ for all $t\to +\infty$, namely,

We emphasize that the Łojasiewicz–Simon type inequality played an essential role in obtaining the required $L^1$-integrability of $\|\partial_t\phi(t)\|_{(H^1(\Omega))'}$ on the unbounded interval $\mathbb{R}^+$ from its $L^2$-integrability on $\mathbb{R}^+$, where the latter is an easy consequence from the energy dissipation (2.11). To achieve such a goal (not limited to the Cahn–Hilliard equation but also for other gradient-like systems), different arguments have been developed in the literature, we may refer to, for instance, ^{[77,94,120,121]} and the references cited therein.

Later in ^[65], the authors studied problem (1.1)–(1.4) with $M = 1$ and a logarithmic type potential like (2.2). They derived an extended Łojasiewicz–Simon type inequality (see [65,Proposition 6.3]) and applied it to show that every global solution converges to a single equilibrium as $t\to+\infty$ in $H^{2r}(\Omega)$ for $r\in (0, 1)$. In this case, we note that the regularity of stationary solutions, in particular, the property of strict separation from pure states (which is indeed a direct consequence from the observation that every solution of the stationary Cahn–Hilliard equation solves the viscous Cahn–Hilliard equation), plays a crucial role in the derivation of the Łojasiewicz–Simon inequality involving a singular potential. In a similar sprit, the authors of ^[97] presented a direct constructive descent method for finding minimizers of nonconvex functionals via the Łojasiewicz–Simon approach and applied it to phase separation problems in multicomponent systems as well as image segmentation. Recently, in ^[55] the authors extended the convergence result for global weak solutions to the Cahn–Hilliard–Oono equation (2.24) with a logarithmic potential subject to (1.2)–(1.4), see [55,Theorem 7.1]. Besides, we would like to mention that the above convergence results for problem (1.1)–(1.4) can be extended to the case with non-constant mobilities, provided that $M(\phi)$ is non-degenerate.

3.   Cahn–Hilliard equation with dynamic boundary conditions

The influence of boundaries (solid walls) on the phase separation process of binary mixtures has attracted a lot of attentions of scientists. For instance, the occurring structures of a binary polymer mixture during the phase separation process may get frozen by a rapid quench into the glassy state and micro-structures at surfaces on small length scales can be produced ^[95]. In order to describe the short-range interactions between the wall and both components of the mixture, suitable surface free energy functional should be introduced into the system (see e.g., ^[95])

where $\nabla_\Gamma$ stands for the tangential (surface) gradient operator defined on the boundary $\partial\Omega$ and $G$ is a surface potential function that may take different forms in various physical contexts. The coefficient $\kappa\geq 0$ is related to the influence of spatial order parameter fluctuations on $\partial \Omega$ (e.g., the surface diffusion). When $\kappa = 0$, the model is closely related to the evolution of a free interface in contact with the solid boundary, namely, the well-known moving contact line problem ^{[19,20,73,106,107]}.

Like in Section 2, the Cahn–Hilliard equation for phase separation processes with boundary effects still has to be supplemented by two boundary conditions, with natural considerations from the mass conservation and energy dissipation (or in other words, energy balance). However, the crucial difference is that one now has to take into account the non-trivial boundary effects driven by the surface energy (3.1). Different types of dynamic boundary conditions have been proposed and analyzed in the literature. In what follows, we shall review some recent progresses in this direction. Before proceeding, we just remark that dynamic boundary conditions can be found in various physical problems (e.g., heat conduction with heat source on the boundary, kinetic motion on the boundary of a vibrating object etc), see for instance, ^[92,93,115] and the references cited therein. In particular, we refer to ^{[74,76,105,117]} for the Caginalp phase-field system with dynamic boundary conditions accounting for the non-isothermal phase transition process of a two-phase material.

Remark 3.1. In the remaining part of this section, since we do not consider the asymptotic behavior with respect to the parameter $\epsilon$ in $E^{\mathrm{bulk}}$ (recall (2.1)), without loss of generality, we simply set $\epsilon = 1$.

3.1.   Dynamic boundary conditions of Allen–Cahn type

The first set of boundary conditions we would to mention are those proposed in ^[95,123] (see also ^[66] for a derivation from the semi-infinite Ising model with Kawasaki spin exchange dynamics). To ensure the conservation of mass in the bulk (cf. (2.8)), one can again assume the following no-flux boundary condition (for the sake of simplicity, hereafter we set the mobility $M = 1$)

On the other hand, to guarantee that the Cahn–Hilliard system tends to minimize its total free energy, the following variational boundary condition was proposed:

where $\Gamma_s > 0$ denotes the surface kinetic coefficient, $\Delta_\Gamma$ stands for the Laplace–Beltrami operator on $\partial \Omega$ and the possible interaction with the bulk part is presented by the term $\partial_{\bf{n}}\phi$. The above boundary condition is usually referred to as a dynamic one, since it contains the time derivative of the order parameter. On the other hand, it can be viewed as a relaxation dynamics ($L^2$-gradient flow) of the surface free energy $E^\mathrm{surf}$ on $\partial \Omega$. Then one can verify that the following energy dissipation relation holds:

When $\kappa = 0$, the dynamic boundary condition (3.3) was proposed to describe some other physically relevant situation for the dynamics of fluid–fluid free interface at the solid boundary (see e.g., ^[20]) such that the contact angle on $\partial \Omega$ turns out to be time dependent and may deviate from the static one like $\pi/2$ (cf. (1.4)).

3.1.1.   Gradient flow structure

Besides the Lyapunov structure given by (3.4), the Cahn–Hilliard equation (1.1) subject to boundary conditions (3.2)–(3.3) can actually be interpreted as a gradient flow of the total energy $E^{\mathrm{total}}(\phi) = E^{\mathrm{bulk}}(\phi)+E^\mathrm{surf}(\phi)$ in the space $\dot{H}^{-1}(\Omega)\times L^2(\partial\Omega)$ with respect to the following inner product:

Roughly speaking, the gradient flow is of a mixed type such that it is $\dot{H}^{-1}$ in the bulk and $L^2$ on the boundary (see [111,Section 3] for further details).

3.1.2.   Well-posedness

Now we review some presentive works on the well-posedness of the initial boundary value problem

Let us first focus on the case with $\kappa > 0$. The first result on the existence and uniqueness of global strong solutions to problem (3.5) was obtained in ^[134]. There the authors considered the following specific choice of free energies

where the constant $g_s > 0$ accounts for a modification of the effective interaction between the components at the solid wall, and $h_s\neq 0$ describes the possible preferential attraction of one of the two components by the wall ^[123]. To overcome the mathematical difficulties due to the presence of the dynamic boundary condition as well as the Laplace–Beltrami operator on the boundary, they introduced an approximate problem related to the phase-field system of Caginalp type

for some positive parameter $\varepsilon > 0$. Here, we note that the system (3.7) is slightly different from the exact form of the approximate problem in [134,Section 2]. The reason is that we just want to describe the main idea but not to get into some technical details therein. The approximate problem (3.7) can be locally solved by first analyzing its linearized problem and then applying the contraction mapping theorem (see [134,Theorem 3.1]). After deriving some uniform-in-time a priori estimates, the authors were able to extend the unique local solution to be a global one (see [134,Theorem 4.1]). Finally, with the help of uniform a priori estimates that are independent of $\varepsilon$, they passed to the limit as $\varepsilon\to 0^+$, to obtain the existence of global strong solutions to the original problem (3.5).

Uniqueness of solutions can be proved by using the standard energy method, however, it was not clear whether the solution obtain in ^[134] defines a $C_0$-semigroup in the energy space (cf. (2.1), (3.1) and note that here $\kappa > 0$):

To handle this issue, in a subsequent paper ^[132], the authors studied the maximal $L^p$-regularity of the system (see [132,Theorem 2.1]) and provided an alternative proof on the global well-posedness of problem (3.5) (see [132,Theorem 4.1]), which further implies that the solution indeed defines a $C_0$-semigroup in $\mathcal{V}^1$ as expected. Based on this fact, they were able to prove the existence of a global attractor in suitable phase spaces, see [132,Theorem 5.1,Theorem 5.2]. Later, the maximal $L^p$-regularity approach was also used in ^[133] to analyze the non-isothermal Cahn–Hilliard equation with dynamic boundary conditions. For further applications to parabolic problems with boundary dynamics of relaxation type, we refer to ^[90].

A third approach to study problem (3.5) was given in ^[129], where the authors considered the viscous Cahn–Hilliard equation for $\alpha\geq 0$:

where $F, G\in C^3(\mathbb{R})$ are some general (but regular) potentials with arbitrary growth and satisfy the following dissipative conditions:

They introduced a new variable on the boundary, i.e., by taking the trace of the phase function $\phi$ such that

and then treated the dynamic boundary condition as a separate evolution equation (parabolic when $\kappa > 0$) on $\partial \Omega$, namely,

Within this convenient framework, the authors proved the existence and uniqueness of global solutions to problem (3.10) by using the Leray–Schauder principle, see [129,Theorem 2.1,Corollary 2.2,Theorem 2.2]. Moreover, they constructed a robust family of exponential attractors (as $\alpha\to 0^+$) for the semigroups associated with problem (3.10) in suitable phase spaces under the constraint due to mass conservation (see [129,Theorem 3.1]). As a byproduct, the result therein also implies that the global attractor obtained in the previous work ^[132] has finite fractal dimension.

The analysis for problem (3.5) with singular potentials are more involved, since the combination of dynamic boundary conditions with singular potentials can produce additional strong singularities on the corresponding solutions close to the boundary (see ^[130]). In this direction, the first result on the existence and uniqueness of global weak solutions was obtained in ^[113], within the following general setting

where $\alpha, \kappa\geq 0$, see [113,Theorem 1,Theorem 2]. In the above system, the authors introduced a suitable decomposition of the singular functions $F'$ and $G'$ such as

respectively, where $\beta$ and $\beta_\Gamma$ are monotone and possibly non-smooth, while $\pi$ and $\pi_\Gamma$ are some regular (Lipschitz) perturbations. In ^[113], they aimed to keep the form of nonlinearities as general as possible. For instance, $\beta$ and $\beta_\Gamma$ are allowed to be essentially arbitrary, with some compatibility conditions such that $\beta$ grows faster than $\beta_\Gamma$ (i.e., the bulk potential plays a dominating role) and that the other boundary contributions satisfy a specific sign condition (see [113,Remark 6]). Well-posedness results under certain alternative assumptions for $\beta$, $\beta_\Gamma$, such as some growth conditions, were also discussed in [113,Theorem 3]. The proofs for the existence results in ^[113] relied on the Yosida type regularization of the possibly singular potentials, a suitable Faedo–Galerkin scheme (via eigenfunctions of the elliptic problem $-\Delta u = \lambda u$ subject to the homogeneous Neumann boundary condition) together with the compactness method (based on suitable a priori estimates performed on the approximate solutions). This argument indeed provides the fourth approach to handle problem (3.5) in a general setting.

Some further investigation was performed in ^[130]. We recall that the result obtained in ^[113] yields the existence of global weak solutions if the (regular) surface nonlinearity $G'$ has the "right" sign at the singular points of the bulk nonlinearity $F'$ (i.e., at the pure states), that is, $\pm G'(\pm 1) > 0$. However, when the sign condition is violated, it was shown in [130,Remark 6.2] that one may not be able to have classical solutions even in the one dimensional case. To overcome this difficulty, the authors of ^[130] introduced a suitable notion of variational solutions (see [130,Definition 3.1]) and proved its existence as well as uniqueness by using proper approximations of the singular potential (see [130,Theorem 3.2]). Then they proved the existence of global and exponential attractors in [130,Theorem 5.2]. Besides, connections between the variational solution and the solution in the usual sense of distributions were established in ^[130]. The authors also discussed the possible separation of the solutions from the singularities of the singular bulk potential (see [130,Section 4]). They showed that those variational solutions are Hölder continuous in space and will become solutions in the usual sense if they do not reach the pure states on the boundary. This property can be guaranteed if either the sign condition as in ^[113] holds or the singularity of the bulk potential $F$ is strong enough (unfortunately, it is not satisfied by the logarithmic one (2.2)). We refer to ^{[27,41,58,130]} for detailed discussions.

In ^[113,130], the authors mainly treated the case that the bulk potential $F$ plays a dominating role. Motivated by the study for the Allen–Cahn equation with singular potentials and dynamic boundary conditions in ^[69], the authors of ^[85] considered the opposite side of the compatibility condition for potential functions such that the boundary potential $G$ is now the leading one (see [85,(2.11)]). For problem (3.9) with $\alpha \geq 0$ and some general assumptions on the potentials that cover all these typical cases (2.2), (2.3) and (2.4), the authors proved global existence, uniqueness and regularity results on its solutions (see [85,Theorem 2.2,Theorem 2.3,Theorem 2.4,Theorem 2.6]). The argument therein implies that in the case of a dominating boundary potential, the analysis can be somewhat simplified and it indeed allows for a unified treatment for the initial boundary value problem (3.9). As far as some extensions are concerned, we refer to ^[79] for the Cahn–Hilliard system with dynamic boundary conditions and mass constraint on the boundary, and to ^[86] for the study of a nonstandard viscous Cahn–Hilliard system with dynamic boundary condition.

When $\kappa = 0$, i.e., the surface diffusion on $\partial \Omega$ is absent, we recall that existence and uniqueness of global weak solutions to problem (3.5) have already been discussed in ^[113]. We also refer to ^[73] for the analysis of problem (3.5) in a different context, which accounts for the evolution of fluid-fluid free interfaces along the solid boundary in the "slow" dynamics such that the effect of flow can be neglected. The authors proved well-posedness of the system with regular potentials (see [73,Theorem 2,Theorem 3]) and investigated its sharp interface limit via the method of matched asymptotic expansion. Besides, the dynamics of the contact point and the contact angle were described and the theoretic results were compared with numerical simulations.

Then a natural question arises, when the surface diffusion in the dynamic boundary condition is vanishing as $\kappa\to 0^+$, whether the solution to problem (3.5) with surface diffusion will converge to the solution of problem (3.5) without surface diffusion? In ^[82], the authors gave an affirmative answer for a wide class of potentials including (2.2), (2.3) as well as (2.4) (see [82,Theorem 2.2]). Their result indicates that the solution of the limiting problem with $\kappa = 0$ will lose some spatial regularity due to the absence of the surface diffusion on $\partial\Omega$ as expected. For instance, the normal derivative $\partial_{\bf{n}}\phi$ on the boundary in general is no longer a function but an element in a dual space. Only under some specific conditions it (as well as the boundary condition (3.3)) can hold almost everywhere (see [82,Theorem 5.1]).

3.1.3.   Long-time behavior

When the long-time behavior of global solutions to problem (3.5) is concerned, we refer to ^{[27,41,114,129,130,132]} for extensive studies within the theory of global and exponential attractors, see also ^{[70,71,100,101,108,109]} for related results on some extended systems involving thermal or memory effects. The role of surface diffusion in dynamic boundary conditions for parabolic and elliptic equations was also discussed in ^[103] within the attractor theory.

Thanks to the energy dissipation relation (3.4), we see that the dynamics of the system can be viewed as a mixing of the Cahn–Hilliard type in the bulk combined with the Allen–Cahn type on the boundary. In particular, we expect that the global solution of problem (3.5) will converge to a single equilibrium as $t\to+\infty$ under suitable assumptions on the potentials $F$ and $G$. When $F$, $G$ take the specific form in (3.6), in ^[144], the authors proved the convergence result (see [144,Theorem 1.1]) by applying an extended Łojasiewicz–Simon type inequality that involves a boundary term (see [144,Lemma 3.4]). In ^[78], the authors obtained the convergence of global solutions to problem (3.5) with a general regular potential $F$ being analytic in the phase function and under suitable growth assumptions (see [78,Theorem 2.3]). They derived a Łojasiewicz–Simon type inequality in the dual space (see [78,Proposition 6.6]) by adapting the abstract result given in [77,Corollary 3.11]. Moreover, they proved the convergence to a single equilibrium of global solutions to the phase-field model of Caginalp type with dynamic boundary conditions (cf. (3.7)). In this aspect, we also refer to ^[143] for the convergence result of an extended phase-field system with hyperbolic relaxation but without boundary diffusion. More general case was treated in ^[114], where the authors considered a general singular bulk potential (including (2.2)) together with a regular boundary potential, in the framework of ^[113]. They extended the argument in ^[78] and derived an extended Łojasiewicz–Simon inequality that also allows nonlinear term on the boundary, see [114,Proposition 6.1] (recall that in ^[78,144] only linear boundary conditions were considered). Then they proved the convergence to a single equilibrium for global solutions as $t\to+\infty$. For related results on some extended Cahn–Hilliard systems, we may refer to ^[70,71,106] and the references cited therein.

We also mention ^[104] for the study on a nonlocal version of the Cahn–Hilliard equation characterized by a fractional diffusion operator which is subject to a fractional dynamic boundary conditions. For the case with regular potentials, the author proved global well-posedness, regularity of solutions (see [104,Theorems 3.3,3.4,3.5]) and the existence of an exponential attractor (see [104,Theorem 4.3]), which yields the existence of a global attractor with finite fractal dimension (see [104,Corollary 1]).

3.2.   Dynamic boundary conditions of Cahn–Hilliard type (Ⅰ)

For simplicity, hereafter we consider the Cahn–Hilliard equation (1.1) with $M = 1$. In ^[98], the author considered the phase separation in a binary mixture confined to a bounded domain with porous walls, which may be (semi) permeable. In this case, due to the possible mass transfer on and through the boundary $\partial\Omega$, the homogeneous Neumann boundary condition (3.2) for the bulk chemical potential $\mu$ should be modified in a proper way. Alternatively, the author of ^[98] proposed the following Wentzell type boundary condition for $\mu$, that is

with some coefficients $b > 0, c\geq 0$ (cf. ^[93,115] for the Wentzell boundary condition for the heat equation and wave equation). Then for sufficiently smooth solutions $(\phi, \mu)$, (3.12) is equivalent to the following dynamic boundary condition

Under the above choice, we no longer have mass conservation in the bulk (unless letting $b\to+\infty$), but a new conservation law for the "total mass" defined in the measure space $(\overline{\Omega}, {\mathrm{d}}\nu) = (\Omega, {\mathrm{d}} x)\oplus(\partial\Omega, {\mathrm{d}} S/b)$ (see ^[98,115]):

The boundary condition (3.12) describes the situation that there is certain mass source (or leak) on the boundary and thus the system undergoes a different diffusive process. When $c = 0$, we see that the total mass is conserved along the evolution such that

Namely, the wall is non-permeable, but certain mass exchange between the bulk and the boundary is allowed. On the other hand, in order to guarantee that the system tends to minimize its total free energy, the author proposed the following variational boundary condition in ^[98] (with $G'$ being a linear function of $\phi$):

As a consequence, it holds

Later in ^[116], the authors considered the phase separation process in a bounded domain with non-permeable walls in some more general setting. Requiring that the total free energy is decreasing in time and the total mass is conserved (see (3.14)), they derived the following set of boundary conditions via a variational principle (see also ^[75]),

for some $b > 0$, $\sigma\geq 0$ and $\kappa\geq 0$. Here, the term $\Delta_\Gamma \mu$ in (3.17) corresponds to some regularizing effect for the chemical potential $\mu$ on the boundary. For $\sigma > 0$ and $\kappa > 0$, (3.17) is often referred to as the Cahn–Hilliard type dynamic boundary condition in the literature. When $\sigma = 0$, (3.17) simply reduces to the Wentzell boundary condition (3.13) with $c = 0$. In ^[116], the authors indeed considered a more general situation such that the parameter $b$ is replaced by a positive function $w^{-1}$: $w\in L^\infty(\partial\Omega)$, $0 < w_*\leq w(x)\leq w^*$ for almost everywhere $x\in \partial\Omega$, where $w_*, w^*$ are given constants. With this generalization, they were able to handle the situation when the "boundary mass" is linked to the variable $\phi$ in a different way with respect to the "bulk mass", for instance, the case when the dynamic boundary condition arises as an approximation of a thin diffusive layer occupied by a different material (see [116,Section 1]). From (3.17) and (3.18), it follows that the new mass conservation property (3.14) is satisfied and the following energy dissipation relation holds (cf. (3.16))

3.2.1.   Gradient flow structure

As pointed out in ^[111], the Cahn–Hilliard equation (1.1) (with $M = 1$) subject to (3.17)–(3.18) admits a gradient flow structure. For simplicity of the presentation, we set $b = 1, \sigma > 0$. Let us consider the elliptic problem of Wentzell Laplacian

Applying the Lax–Milgram theorem, one can see that the above problem admits a unique weak solution $u = \mathcal{W} (f)\in \mathcal{V}^1$ with the constraint $\int_\Omega u\, {\mathrm{d}} x+\int_{\partial\Omega} u\, {\mathrm{d}} S = 0$, if the two nonhomogeneous terms on the right-hand side satisfy $f = (f_\Omega, f_\Gamma)\in (\mathcal{V}^1)'$ and $\langle f_\Omega, 1\rangle_{(H^1(\Omega))', H^1(\Omega)}+\langle f_\Gamma, 1\rangle_{(H^1(\partial\Omega))', H^1(\partial\Omega)} = 0$. With the help of the solution operator $\mathcal{W}$, Eq (1.1) subject to (3.17)–(3.18) can be viewed as a gradient flow of the total energy $E^{\rm{total}}$ with respect to the inner product

Similar conclusions could be drawn for more general choices of parameters.

3.2.2.   Well-posedness and long-time behavior

Assume that $\Omega\subset \mathbb{R}^d$, $d = 1, 2, 3$, is a smooth bounded domain. Without loss of generality, let us write the resulting initial boundary value problem in the following form:

Note that in the one dimensional case, the terms involving $\Delta_\Gamma$ do not appear in any of the boundary conditions.

We first review some known results for the case $\sigma = 0$. The first results on existence and uniqueness of global strong solutions to problem (3.20) with a regular bulk potential $F$ and a quadratic boundary potential $G$ was obtained in ^[98] by employing some classical methods such as the contraction mapping principle and the energy method (see [98,Theorem 11,Theorem 12]). The proof therein extends the idea in ^[134] by investigating an approximate system of Caginalp type. Besides, a comparison between two solutions of the systems (3.5) and (3.20) was made in [98,Theorem 14], which provides an estimate on the difference between the two solutions on certain finite time interval, explicitly in terms of the parameters $b$ and $\Gamma_s$. Later in ^[99], the author improved his result by using a different approach, i.e., the Faedo–Galerkin method, to prove the existence and uniqueness of a global solution to the same problem, but under some more general assumptions on the bulk potential $F$, see [99,Theorem 4.1]. The result implies that the solution defines a continuous semigroup on a suitable phase space, which enables him to investigate its long-time behavior by proving the existence of an exponential attractor (and thus of a global attractor with finite dimension), see [99,Theorem 5.1,Theorem 5.6]. In ^[102], the author further studied the asymptotic behavior as $b\to+\infty$ and constructed a robust family of exponential attractors for the problem under the same assumptions on the potential as in ^[99].

Convergence of global solutions to a single equilibrium as $t \to +\infty$ for problem (3.20) with $\sigma = 0$ was first analyzed in ^[142]. Under the assumption that $c, \kappa > 0$ and $F$ is analytic with respect to the phase function, the author proved the expected convergence result [142,Theorem 1.1] by means of an extended Łojasiewicz–Simon type inequality with boundary term (see [142,Lemma 3.5]). Higher order estimates on the convergence rate were also obtained in ^[142] by combining the Łojasiewicz–Simon approach (cf. ^[119]) together with the energy method. The mass conserved case $c = 0$ is somewhat more involved, since the boundary contribution in the energy dissipation of (3.16) vanishes. To overcome the difficulty due to a weaker dissipation that only comes from the bulk, in ^[110] the authors applied a generalized Poincaré inequality to derive an improved Łojasiewicz–Simon type inequality subject to certain mass constraint (see [110,Lemma 4.1]), and then achieved the convergence result [110,Theorem 2.4]. Moreover, they were able to prove the convergence of global solutions to steady states in the limiting case $b = \infty$.

Next, the case $\sigma > 0$, $c = 0$ was studied in ^[116] in a rather general setting, i.e., the bulk potential $F$ can be a singular one in a wide class and the constant $b$ can be a bounded positive function. There the authors overcame difficulties similar to those encountered in ^[130] and proved existence, uniqueness as well as regularity of global weak solutions (defined in a suitable sense) and studied their long-time behavior, including the existence of a compact global attractor and the convergence to a single equilibrium as $t\to +\infty$. The arguments in ^[116] were rather involved, since the authors tried to deal with the nonlinear terms under general assumptions, see [116,Section 3] for detailed discussions. Concerning a different approach, we refer to ^[122] for existence and uniqueness of solutions to problem (3.20) with regular potentials on a domain with either permeable or non-permeable walls, which were obtained by applying the general theory of maximal $L_p$ regularity of relaxation type in ^[90].

We note that the boundary nonlinearity $G$ in ^[116] was assumed to be a regular function. The case with a possibly singular and dominating boundary potential was considered in ^[80]. The authors introduced the unknowns on the boundary

and considered the following alternative formulation of (3.20) (with some specific choice of parameters)

where $f, f_\Gamma$ are some given functions. Working with a general setting of bulk and boundary potentials including the typical types (2.2), (2.3) and (2.4), they proved existence and uniqueness (indicated by a continuous dependence estimate) of global weak solutions, see [80,Theorem 2.1,Theorem 2.2]. Refined regularity for weak solutions and existence of strong solutions were obtained as well, see [80,Theorem 4.2]. The proof of the existence result relies on the Yosida regularization of maximal monotone graphs and the introduction of viscous terms in the equations for both bulk and boundary chemical potentials $\mu$, $\mu_\Gamma$ in (3.21). Then the solvability of the approximate problem can be deduced from the abstract theory of doubly nonlinear evolution inclusions ^[89]. After obtaining uniform estimates for the approximate solutions, the authors of ^[80] were able to pass to the limit and conclude the existence of global weak solutions. See ^[81] for an extension to the Cahn–Hilliard equation on the boundary with bulk condition of Allen–Cahn type. We also refer to ^[87,88] for extended results on well-posedness and long-time behavior for the viscous regularization of problem (3.20) with an additional convection term. The proof for the existence of solutions in ^[87] was based on a similar regularization of maximal monotone graphs but a different way of approximation for the equation via the Faedo–Galerkin method.

In the recent work ^[96], for problem (3.20) with $\sigma > 0$, $c = 0$ and physically relevant singular (e.g., logarithmic) potentials, global regularity of weak solutions was established, see [96,Theorem 2.1]. When the spatial dimension is two, the authors showed the instantaneous strict separation property such that for arbitrary positive time any weak solution stays away from the pure phases $\pm 1$, while in the three dimensional case, although the instantaneous separation property remains open, an eventual separation property for large time was derived. As a consequence, they proved that every global weak solution converges to a single equilibrium as $t\to+\infty$ (see [96,Theorem 2.2]), by applying a suitable Łojasiewicz–Simon type inequality. We expect that a similar result could be obtained for the case $\sigma = 0$, $c\geq 0$.

3.3.   Dynamic boundary conditions of Cahn–Hilliard type (Ⅱ)

From the previous discussions, we see that under all choices of boundary conditions ((3.2) with (3.3), (3.13) with (3.15), or (3.17) with (3.18)), the Cahn–Hilliard equation (1.1) satisfies two basic physical constraints, that is, the mass conservation and the energy dissipation. Among them, the boundary conditions (3.2), (3.13) and (3.17) are proposed to keep suitable mass conservation property in the physical domain (see (2.8), (3.14)), while the so-called variational boundary conditions (3.3), (3.15), (3.18) are chosen in a phenomenological way so that the validity of some specific energy dissipation relation is guaranteed (see (3.4), (3.16), or (3.19)). We note that (3.3), (3.15) and (3.18) can be viewed as some sufficient conditions for the energy dissipation of the evolution system, however, such choices may not be unique.

3.3.1.   Model derivation

In ^[37], the authors derived a different type of dynamic boundary condition for the Cahn–Hilliard equation (1.1), that is,

The derivation is based on an energetic variational approach that combines the least action principle and Onsager's principle of maximum energy dissipation, from which we see that the equation (1.1) subject to (3.22) naturally fulfills three basic physical properties: (1) kinematics: conservation of mass both in the bulk $\Omega$ and on the boundary $\partial \Omega$; (2) energetics: dissipation of the total free energy; (3) force balance: both in the bulk $\Omega$ and on the boundary $\partial\Omega$.

Below we sketch the derivation of (3.22), and refer to ^[37] for further details in a more general setting. In the bulk $\Omega$, $\phi$ is assumed to be a locally conserved quantity that satisfies the continuity equation

where ${\bf{u}}: \Omega \to \mathbb{R}^d$ stands for the microscopic effective velocity (e.g., due to diffusion process etc). We assume that ${\bf{u}}$ satisfies the no-flux boundary condition

Next, we consider some nontrivial boundary dynamics that is assumed to satisfy a local mass conservation law analogous to (3.23) such that

where ${\bf{v}}: \partial\Omega \to \mathbb{R}^d$ denotes the microscopic effective tangential velocity field on the boundary. We note that there is no need to impose any boundary condition on ${\bf{v}}$, since here the boundary $\partial \Omega$ is assumed to be a closed manifold.

Next, for an isothermal closed system, evolution of the binary mixtures is assumed to satisfy the following energy dissipation law

where

For instance, $E^{{\mathrm{bulk}}}$ and $E^{\mathrm{surf}}$ are given by (2.1) and (3.1), respectively. On the other hand, the rate of energy dissipation $\mathcal{D}^{\mathrm{total}}$ is chosen as

which also consists of contributions from the bulk and the boundary. Here, we assume that

where $M_{\mathrm{b}}$, $M_{\mathrm{s}}$ are some positive mobility functions.

Finally, in order to derive a closed system of partial differential equations, it remains to determine the microscopic velocities ${\bf{u}}$, ${\bf{v}}$ in equations (3.23) and (3.25). Using the principle of least action and Onsager's maximum dissipation principle, we are able to derive the conservative and dissipative forces according to the free energy (3.27) and the dissipation functional (3.28). Then by the force balance relation (Newton's second law), we obtain

where the bulk and boundary chemical potentials are given by

Solving ${\bf{u}}$, ${\bf{v}}$ from (3.30) and inserting them back into (3.23) and (3.25), respectively, we arrive at the following Cahn–Hilliard system subject to a new class of dynamic boundary condition:

When the mobilities are positive constants, for instance, $M_{{\mathrm{b}}} = M_{{\mathrm{s}}} = 1$, it was shown in ^[111] that problem (3.31) can be regarded as a $H^{-1}$ type gradient flow equation of the total free energy $E^{\mathrm{total}}$ both in the bulk and on the boundary, with respect to a suitable inner product on certain function space (see [111,Section 3] for details).

3.3.2.   Well-posedness

Inspired by ^[129], it will be convenient to view the trace of the order parameter $\phi$ as an unknown function on the boundary. After introducing the new variable

the initial boundary value problem of the Cahn–Hilliard system (3.31) can be written in the following form (taking $M_{\rm{b}} = M_{\rm{s}} = 1$ for simplicity):

Assume that $\Omega\subset \mathbb{R}^d$, $d = 2, 3$, is a smooth bounded domain. Well-posedness of problem (3.32) was first analyzed in ^[37] when $F$ and $G$ are suitable regular potentials (including the typical choice like (2.3) for $F$ and the physically relevant surface potential in ^[20] for $G$). The one dimensional case would be easier since the boundary dynamics is trivial such that the terms involving the operator $\Delta_\Gamma$ simply vanish. The proof therein was inspired by the argument in ^[129]. First, for the problem with surface diffusion (i.e., $\kappa > 0$), the authors introduced a regularization by adding viscous terms in both of the bulk and boundary chemical potentials. This leads to a viscous Cahn–Hilliard equation subject to a dynamic boundary condition of viscous Cahn–Hilliard type, which can be solved by using the contraction mapping principle (see [37,Proposition 4.1]). After deriving global-in-time a priori estimates that are independent of the approximating parameter as well, they obtained the existence of global weak (and strong) solutions by finding a convergent subsequence of the approximate solutions after passing to the limit, see [37,Theorem 3.1].

For the case without surface diffusion, i.e., $\kappa = 0$, the existence result can be achieved by deriving uniform estimates independent of the parameter $\kappa$ and then taking the limit as $\kappa \to 0^+$ (i.e., the vanishing surface diffusion limit), see [37,Theorem 3.2]. In both cases, the uniqueness of solutions can be proved by a standard energy method. It is worth mentioning that the solution will lose some spatial regularity due to the absence of the surface diffusion (cf. ^[82] for a similar situation for problem (3.5)). In order to guarantee that the boundary equation $\mu_\Gamma = -\partial_{\bf{n}}\phi + G'(\psi)$ is satisfied for weak solutions in the usual sense, that is, almost everywhere on $\partial\Omega \times (0, T)$, an additional geometric assumption was imposed in [37,Theorem 3.2]. This restriction was later removed in ^[111], where the authors introduced a weaker notion of the "weak solution" such that the bulk and boundary chemical potentials satisfy the following weak form (see [111,Definition 4.1])

for any test function $\eta\in L^2(0, T; H^1(\Omega))\cap L^\infty(\Omega\times (0, T))$ with $\eta|_{\partial\Omega}\in L^\infty(\partial\Omega\times(0, T))$. Namely, the troublesome boundary term $\partial_{\bf{n}} \phi$ does not appear in an explicit way. Then they proved existence and uniqueness of global weak solutions by an implicit time discretization based on the gradient flow structure of problem (3.32), see [111,Theorem 4.3]. Unfortunately, their argument does not work for singular potentials like the logarithmic potential (2.2) or the obstacle potential (2.4).

Well-posedness of problem (3.32) with singular potentials was first established in ^[84]. For the case $\kappa > 0$, the authors treated the initial boundary value problem in a wide class of nonlinearities (covering (2.2), (2.3) and (2.4)) with the compatibility condition that the boundary potential plays a dominating role. They proved existence and uniqueness of global weak solutions (see [84,Theorem 2.3,Theorem 2.4]) as well as existence of a unique global strong solution (see [84,Theorem 4.1]). The proofs therein relied on several approximations of the original problem such as the introduction of viscous regularizations in the chemical potentials, the Yosida approximation for maximal monotone graphs and an implicit time discretization scheme for the bulk-boundary coupled system. They first proved the existence of a discrete solution by taking advantage of the general maximal monotone theory. After deriving a number of uniform estimates, they were able to conclude the existence results by performing limiting procedures with respect to the time step first and then the approximating parameters.

3.3.3.   Long-time behavior

Concerning the long-time behavior of problem (3.32), the authors made a preliminary investigation for the case $\kappa > 0$ in ^[128] within the framework of infinite dimensional dynamical systems. For the system with regular potentials, they proved existence of exponential attractors, which also yields the existence of a global attractor with finite fractal dimension (see [128,Theorem 2.1]), while for the system with singular potentials, they showed the existence of a global attractor in a suitable complete metric space (see [128,Theorem 2.2]). The result is less satisfactory for the latter, because of the possible singularity of the bulk and boundary potentials at the pure phases $\pm 1$ and its interplay with the dynamic boundary condition (cf. ^[130]). To overcome this difficulty, one crucial step is to establish the strict separation property of the solution (cf. ^[96] for problem (3.20)). Besides, the case without boundary diffusion ($\kappa = 0$) turns out to be more difficult, since it is not clear so far whether the solutions to problem (3.32) can define a $C_0$-semigroup in a suitable phase space.

The long-time behavior of a single global solution was first analyzed in ^[37]. With the additional assumption that the regular potentials $F$ and $G$ are real analytic functions on $\mathbb{R}$, the authors proved that every global bounded weak/strong solution to problem (3.32) will converge to a single equilibrium as $t\to+\infty$ and provided an estimate on the convergence rate, see [37,Theorem 3.3]. The conclusion was again achieved by employing an extended Łojasiewicz–Simon type inequality [37,Lemma 6.3]. We note that the above convergence results are valid for both cases with or without boundary diffusion, i.e., for $\kappa\geq 0$. Moreover, in presence of the boundary diffusion, by applying the Łojasiewicz–Simon approach in a different way, the authors were able to give a further characterization on the Lyapunov stability of steady states (e.g., local energy minimizers) that may be non-isolated, see [37,Theorem 3.4]. Whether the corresponding convergence or stability results hold for solutions to the system (3.32) with singular potentials remains an open question, since it essentially relies on the regularity of solutions, in particular, the strict separation property from the pure states $\pm 1$.

At last, we remark that all the available results for problem (3.32) were obtained for the case with constant mobilities. It will be interesting to study the case with non-constant mobilities that are concentration dependent or even degenerate.

3.4.   Some extended models

In the last part, we mention some extensions of the Cahn–Hilliard equation with dynamic boundary conditions that have been studied in the recent literature.

We note that in the formulations (3.10), (3.21) and (3.32), some "strong" relations were imposed on the phase-field function $\phi$ or the chemical potential $\mu$ via the trace operator, namely, in terms of some nonhomogeneous Dirichlet boundary conditions. To provide a more general description of the interactions between the materials in the bulk and the materials on the boundary, extended models with certain relaxed coupling relations between the bulk and boundary variables (e.g., via the Robin type boundary conditions) were introduced and analyzed in ^{[112,124,125]}.

In ^[124], the authors considered the following bulk-boundary coupling system that can be regarded as an extension of (3.32):

for some constants $\alpha\neq 0$, $\beta\in \mathbb{R}$. The affine linear transmission condition between the bulk and boundary phase-field variables (i.e., $\phi = \alpha \psi + \beta$ on $\partial\Omega$) was introduced in order to account for some non-trivial boundary interactions. This condition can be further relaxed by imposing a Robin type boundary condition such that

for some $K > 0$ and a function $H\in C^2(\mathbb{R})$. Under suitable assumptions, the authors proved existence and uniqueness of global weak solutions to problems (3.33) and (3.34) (see [124,Theorem 2.1,Theorem 2.2]). Besides, for the special case $H(s) = \alpha s+\beta$, they showed the weak convergence of solutions as the parameter $K\to 0^+$, and derived an error estimate between solutions of the two models, see [124,Theorem 2.3]. Analysis on some similar extended models for the Allen–Cahn equation with dynamic boundary conditions can be found in ^[83,126].

An analogous Robin type relaxation can be imposed on the chemical potentials, which account for possible reactions between the materials in $\Omega$ and on $\partial \Omega$. In ^[125], the authors studied the following problem

for some constants $\beta\neq 0$ and $L > 0$. The Robin type boundary condition indicates that the mass flux $\partial_{\bf{n}}\mu$ is driven by differences in the chemical potentials and the constant $L^{-1}$ can be interpreted as the reaction rate. Taking $\beta = 1$, this boundary condition also allows people to establish a connection between models (3.21) and (3.32), via the formal limits $L\to 0^+$ and $L\to +\infty$, which correspond to the limit cases of an instantaneous reaction ($L\to 0^+$) and a vanishing reaction rate ($L\to +\infty$), see ^[125] for detailed discussions. Thus, problem (3.35) can be interpreted as an interpolation between (3.21) and (3.32) where the parameter $L$ corresponds to a positive but finite kinetic rate. The authors proved the existence, uniqueness and regularity of weak solutions to problem (3.35) (see [125,Theorem 3.1]) and investigated the asymptotic limits as $L\to 0^+$ and $L\to +\infty$, establishing also convergence rates for these limits, see [125,Theorem 4.1,Theorem 4.2]. Concerning the long-time behavior of global solutions to problem (3.35), we refer to ^[112], where the authors proved the existence of a global attractor and the convergence of global weak solutions to a single equilibrium as $t\to+\infty$ by means of a Łojasiewicz–Simon type inequality, see [112,Theorem 4.9,Theorem 4.13]. Besides, they proved that the global attractor of problem (3.21) is stable with respect to perturbations of the kinetic rate and constructed a robust family exponential attractors for $L\in [0, 1]$ (see [112,Theorem 6.5]).

We note that the results obtained in ^{[112,124,125]} for the extended models are only valid for regular potentials including (2.3). However, singular potentials like the logarithmic potential (2.2) or the obstacle potential (2.4) are not admissible for the arguments therein. This will be an interesting field for the future research.

Acknowledgments

The work of the author was partially supported by NNSFC 12071084.

Conflict of interest

The authors declare there is no conflicts of interest.