Barycentric interpolation collocation method based on Crank-Nicolson scheme for the Allen-Cahn equation

: This paper proposes a numerical scheme for the Allen-Cahn equation that represents a phenomenological model for anti-phase domain coarsening in a binary mixture. In order to obtain a high order discretization in space, we adopt the barycentric interpolation collocation method. The semi-discretized scheme in space is shown to be consistent. The second-order Crank-Nicolson scheme is used for temporal discretization and the simple iteration method is adopted for nonlinear term. Corresponding algebraic system is derived. Numerical examples are demonstrated to validate the e ﬃ ciency of the proposed method.


Introduction
In this paper, we consider the following Allen-Cahn equation with the Dirichlet boundary condition u(x, t) = φ(x, t), x ∈ ∂Ω, t ∈ (0, T ]. (1.2) where Ω is a bounded domain in R d (d = 1, 2). The positive parameter ε representing generally the interfacial width is a very small constant and the reaction term f (u) = F (u) with F(u) = 1 4 (u 2 − 1) 2 .
The function u(x, t) can be viewed as the difference between the concentrations of the two components of the mixture or be defined as the order parameter that indicates the local state of the whole system. For instance, u = 1 and u = −1 stand for two different phases. It is well known that the Allen-Cahn equation possesses energy-decay property associating with the following total free energy functional In fact, taking the inner product for equation (1.1) with u t , we can obtain (u t , u t ) + 1 ε 2 (u 3 − u, u t ) + (∇u, ∇u t ) = 0.
The derivative of the energy E(u(t)) with respect to t gives dE(u(t)) dt = Ω 1 ε 2 (u 3 − u)u t + ∇u∇u t dx = −(u t , u t ) ≤ 0, which implys the total energy is decreasing in time, namely E(u(t 2 )) ≤ E(u(t 1 )), ∀t 1 < t 2 ∈ (0, T ]. (1.4) The Allen-Cahn equation introduced by Allen and Cahn [1] in 1979 is a kind of non-homogeneous semi-linear poisson equation, which often occurs in convection diffusion equations in computational fluid dynamics or reaction-diffusion problems in material science. Hitherto, the Allen-Cahn equation has been widely used to model mathematically in crystal growth [2], image analysis [3] and average curvature-flow rate [4]. In addition, Allen-Cahn equation can also be used to describe competition of biological populations and phenomenon of exclusion [5]. It has become a basic model to study interfacial dynamics and phase transitions in material science [6]. However, the analytical solutions of the phase-field models can not be obtained so that it is greatly necessary and significant to develope stable, efficient and highly accurate numerical methods.
During the past two decades, enormous works have been devoted to numerical simulations and numerical analysis for the Allen-Cahn equation such as finite difference [7][8][9], finite element [10,11], spectral method [12]. Zhai et al. in [8] proposed a novel linearized high-order compact difference method for Allen-Cahn equation with different boundary conditions. Hou et al. presented numerical analysis of fully discretized Crank-Nicolson scheme for fractional-in-space Allen-Cahn equations and proved that the second order temporal discretization scheme preserves the discrete maximum principle in [9]. Feng and Li [10] developed two fully discrete interior penalty discontinuous Galerkin scheme for Allen-Cahn equation. Li et al. [11] introduced an unconditionally energy stable finite element method for Allen-Cahn equation. But the scheme is only second-order accurate in both space and time. Weng and Tang [12] presented two unconditionally stable second-order operator splitting approaches based on Fourier spectral for solving the Allen-Cahn equation. In addition, the operator splitting method combined with the finite difference method or finite element method for Allen-Cahn equation have also been considered in [13,14].
Most of above-mentioned works depend on the mesh generation to solve differential equations, which causes some difficulties especially for high dimension problems with irregular domains. In recent years, meshless methods have captured attention of scholars due to it's advantage of meshless and ability of treating irregular domains. As a novel meshless method , barycentric interpolation method is a good choice for dealing with polynomial interpolations since it is accurate highly, fast, stable and easy on program implementation. Barycentric interpolation method include barycentric Lagrange interpolation and barycentric rational interpolation, both of which can be written the unified formula of barycentric interpolation. The weight functions become extremely big when the nodes are of uniform distribution, which leads to the Runge phenomenon and ruins the merits of the barycentric Lagrange interpolation. Luckily, if the nodes obey the density proportion (1 − x 2 ) − 1 2 such as the families of Chebyshev points that is the simplest clustered point sets, the barycentric Lagrange interpolation has a good numerical stability. Barycentric interpolation also can effecively avoid the accumulated errors caused by difference scheme. Readers can refer to [15,16] for detailed introduction and derivation of the barycentric interpolation. The collocation method based on barycentric interpolation was lately extended to solve various integral equations and partial differential equations including high-dimensional Fredholm integral equation of the second kind [17], Volterra integral equations with weakly singular kernels [18], nonlinear parabolic equations [19] and 2D viscoelastic wave equation [20].
To our knowledge, the works of the error analysis for barycentric interpolation collocation method are comparatively sparse. Recently, Yi and Yao in [21] put forward a steady barycentric Lagrange interpolation method and presented error analysis of system for solving the time-fractional telegraph equation. Based on the above works, we focus on a fully discrete scheme for the Allen-Cahn equation, which is a second-order Crank-Nicolson scheme in time combined with the barycentric interpolation collocation method in space. Moreover, we will give consistency analysis of the semi-discretized scheme in space.
The remainder of the paper is structured as follows: In section 2, the barycentric interpolation collocation method is introduced. Semi-discretized scheme in space and corresponding consistency analysis are carried out in detail for the Allen-Cahn equation in section 3. In section 4, fully discrete scheme based on Crank-Nicolson scheme is presented. Nmerical examples test the accuracy and efficiency of proposed algorithm in section 5 and some conclusions are given in section 6.

Barycentric Lagrange interpolation with the Chebyshev points
Suppose n + 1 distinct interpolation nodes x j be given, together with corresponding a set of real numbers y j ( j = 0, 1, . . . , n). Let p(x) denotes the polynomial of degree at most n, satisfying p(x j ) = y j ( j = 0, 1, . . . , n). As we all known, such polynomial p(x) is unique and can be written in Lagrange form as where L j (x) is the basis function in Lagrange interpolation, which possesses the following properties Defining the barycentric weights by the L j (x) can be expressed as Inserting Eq (2.6) into Eq (2.1), we get Combining Eq (2.3) with Eq (2.6), we have the following identical equation Dividing Eq (2.8) by Eq (2.7) and cancelling the common factor l(x), the barycentric Lagrange interpolation formula for p(x) can be obtained In the paper, we choose the Chebyshev points given by which will ensure that above polynomial interpolation has good numerical stability.

Barycentric rational interpolation
The rational function interpolation based on the idea of mixed function can effectively overcome the instability of interpolation.

Differential matrix
The derivative of p(x) defined as Eq (2.9) with respect to x as (2.18) The first and second order differentiation matrices can be obtained by the following formula [15] (2.20) By mathematical induction, we have the following v-order differential matrix 3. Semi-discretized system and consistency analysis

Semi-discretized scheme based on barycentric interpolation collocation method
Denoting u(x i , y, t) = u i (y, t), i = 0, 1, . . . , M, and fixing variable y, the unknown function u(x, y, t) can be written in barycentric interpolation form where γ l (x) is the basis function in barycentric interpolation on the direction of x. Substituting Eq (3.1) into Eq (1.1) and making Eq (1.1) be identical at the nodes . . . . . .
Denote u i (y j , t) = u i j (t). Similarly, u i (y, t) can be written in barycentric interpolation form where β q (y) is the basis function in barycentric interpolation on the direction of y. Inserting Eq (3.4) into Eq (3.3) and making Eq (3.3) be identical at the nodes y j , j = 0, 1, . . . , N , we get the following ODE system . . . . . .
where β q (y j ) = d 2 β q (y j ) dy 2 = D (2) jq . Introduce the following notations where the sign ⊗ represents the Kronecker product of the matrix, C (2) and D (2) stand for second order differential matrix on nodes x 0 , x 1 , . . . , x M and on nodes y 0 , y 1 , . . . , y N , severally. I M and I N being the identity matrix of M + 1, N + 1 order, respectively.

Consistency analysis for the semi-discretized scheme
In this part, we present consistency estimates of the semi-discretized system (3.6) with the collocation method. Let p(x) is the Lagrange interpolation function approximating u(x). According to interpolation remainder theorem, we have Estimates of (3.7) are considered by the following Lemma that have been proved in Yi [21]. Lemma 1. (see Yi [21]) If u(x) ∈ C m+1 ([a, b]), then the following estimates for functional e(x) defined as (3.7) hold where C 1 , C * 1 and C * * 1 are constant independent of x, e is natural logarithm and l x represents the length of the interval [a, b].
Let p(x, y) is Lagrange interpolation function of u(x, y) satisfying p(x m , y n ) = u(x m , y n ). Defining the error function as we have the following results based on the Lemma 1.
, wherem = max{m, n}. Then the following estimates for functional e(x, y) defined as (3.9) hold where l y represents the length of the interval [c, d]. Similarly, we get Set u(x m , y n , t) be the numerical solution of u(x, y, t), then we have Πu(x m , y n , t) = 0 (3.12) and lim m,n→∞ Πu(x m , y n , t) = 0, (3.13) where Π is differential operator.
For R 1 , we have = e t (x m , y, t) + e t (x m , y n , t).

(3.16)
Applying Theorem 1, we obtain Analogously, we estimate R 2 and R 3 as

18)
|R 3 | = e yy (x m , y, t) + e yy (x m , y n , t) (3.19) For the term R 4 , we have (3.20) Substituting (3.17)-(3.20) into (3.15), this completes the proof. From Theorem 2, we can find that the order of the difference operator Π determines the consistency rate of the semi-discrete scheme.

Fully discretized scheme based on Crank-Nicolson scheme
In the section, we will solve the ODE system (3.6) by the Crank-Nicolson scheme. Partition the interval (0,T] into a uniform mesh with the time step τ = T l : 0 = t 0 < t 1 < . . . < t l = T . Let U k = U(t k ), k = 0, 1, . . . , l. the nonlinear term of (3.6) is expanded using Taylor formula at the node vector U k (k = 0, 1, . . . , l), we get Crank-Nicolson scheme is used for time discretization, we have Inserting U k+ 1 2 = 1 2 U k+1 + U k into Eq (4.2), we obtain fully discretized scheme as (4.3)

Numerical experiments
In this section, we will perform several numerical examples to test the accuracy and energy decline property of the proposed scheme. we consider convergence tests in the first example and three 2D problems with different initial conditions in the next examples.
For convenience, introducing the following error notations where u h and u e denote the numerical solution and the exact solution of the problem, respectively, · ∞ is the L ∞ norm.

Convergence tests
We test numerical convergence order in time for 1D Allen-Cahn equation at first, the initial and boundary condition are given the following analytical solution Fix M = 30, and time step τ is varied. The errors in maximum norm and convergence rates are displayed in Table 1, from which we see that the scheme is second order in time for 1D Allen-Cahn equation. Next, we verify the accuracy and convergence rate in time of the algorithm when applying to the 2D Allen-Cahn equation with the inhomogeneous term g(x, y, t). We consider the following problem on Ω = [−1, 1] 2 × (0, T ] as  Table 2. One can see that both barycentric interpolation collocation methods have high precision, as they only use 8 × 8 mesh nodes can achieve the accuracy of second-order central difference approach with 40 × 40 mesh nodes. In addition, it is also observed that the accuracy of barycentric Lagrange interpolation is slightly higher than that of barycentric rational interpolation especially taking fewer mesh nodes from Table 2. In addition, exponential convergence property of two methods also can be observed from Figure 1.  Choose M = N = 20 and vary the temporal step τ. Table 3 indicates the derived scheme is indeed of second order in time for 2D Allen-Cahn equation. After that, we will show the energy decay property of the 2D Allen-Cahn equation. Define the discrete energy function as where u k i j and E h (u k ) represent the numerical solution and the energy value at the k-th time separately, h i j denotes spatial step size between adjacent Chebyshev nodes. The space diagram of 2D Allen-Cahn equation is shown in Figure 2. And from Figure 3, it is obvious that the energy is decreasing with time.
From Figure 4, we can see clearly that the connected interface splits into two curves at first. Then, the two components of the interface develop circular shapes. Eventually, the diameters of the two particles decrease to zero until they collapse.
(a) t = 0s (b) t = 0.5s (c) t = 2s (d) t = 3s Figure 5. Snapshots of the phase variable u are taken at different time for Example 3.
The coarsening phenomena of Allen-cahn equation with the evolution of time are presented in Figure 6, from which we see that the solution reaches the steady state at t = 0.1s.

Conclusions
In this work, by a combination of the Crank-Nicolson scheme and barycentric interpolation collocation method, an efficient numerical scheme for Allen-Cahn equation is developed. Besides, the consistency analysis for the semi-discretized scheme in space is derived. The numerical examples show that our scheme is second order in time and convergent exponentially in space. Energy dissipation law is also checked numerically. In future work, we plan to give the error estimate of fully discretized scheme and generalize this method to other models.