{'fraction_non_alphanumeric': 0.05857275254865616, 'fraction_numerical': 0.0231490062815364, 'mean_word_length': 4.011456290638868, 'pattern_counts': {'":': 0, '<': 22, '<?xml version=': 0, '>': 18, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 10, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Sharp interface limit in a phase field model of cell motility

We consider a system of two coupled parabolic PDEs introduced in [1] to model motility of eukaryotic cells. We study the asymptotic behavior of solutions in the limit of a small parameter related to the width of the interface in phase field function (sharp interface limit). We formally derive an equation of motion of the interface, which is mean curvature motion with an additional nonlinear term. In a 1D model parabolic problem we rigorously justify the sharp interface limit. To this end, a special form of asymptotic expansion is introduced to reduce analysis to a single nonlinear PDE. Further stability analysis reveals a qualitative change in the behavior of the system for small and large values of the coupling parameter. Using numerical simulations we also show discontinuities of the interface velocity and hysteresis. Also, in the 1D case we establish nontrivial traveling waves when the coupling parameter is large enough.


Introduction
The problem of cell motility has been a classical subject in biology for several centuries. It dates back to the celebrated discovery by van Leeuwenhoek in the 17th century who drastically improved the microscope to the extent that he was able to observe motion of single celled organisms that moved due to contraction and extension. Three centuries later this problem continues to attract the attention of biologists, biophysicists and, more recently, applied mathematicians. A comprehensive review of the mathematical modeling of cell motility can be found in [25]. This work is motivated by the problem of motility (crawling motion) of eukariotic cells on substrates. The network of actin (protein) filaments (which is a part of the cytoskeleton in such cells) plays an important role in cell motility. We are concerned with cell shape dynamics, caused by extension of the front of the cell due to polymerization of the actin filaments and contraction of the back of the cell due to detachment of these filaments. Modeling of this process in full generality is at present a formidable challenge because several important biological ingredients (e.g., regulatory pathways [25]) are not yet well understood.
However, in recent biophysical studies several simplified phase field models have been proposed. Simulations performed for these models demonstrated good agreement with experiments (e.g., [10,13] and references therein). Recall that phase field models are typically used to describe the evolution of an interface between two phases (e.g., solidification or viscous fingering). The key ingredient of such models is an auxiliary scalar field, which takes two different values in domains describing the two phases (e.g., 1 and 0) with a diffuse interface of a small (non zero) width. The advantage of such an approach is that it allows us to consider one set of PDEs in the whole domain occupied by both phases and therefore avoids an issue of coupling two different sets of PDEs in each phase, which is typically quite complicated in both simulations and analysis.
In this work we present rigorous mathematical analysis of the 2D phase field model proposed in [13] that consists of a system of two PDEs for the phase field function and orientation vector with an integral mass conservation constraint. This model can be rewritten in a simplified form suitable for asymptotic analysis, so that all key features of the qualitative behavior are preserved, which can be seen from a comparison of simulations from [13] with our analytical results. First, in [13] the integral mass conservation constraint is introduced in the PDE system by adding a penalization parameter into the double-well potential (formulas (2.2), (2.5)-(2.6) in [13]). We introduce this constraint via a dynamic Lagrange multiplier λ (t) defined below, which provides the same qualitative behavior of solutions. Second, for technical simplicity we drop two terms in the second equation (for polarization) in the phase field system [13], since our analysis shows that these terms can be incorporated with minor changes in both the results and the techniques. Thirdly, in order to study the long term behavior of the system, we perform the diffusive scaling (t → ε 2 t, x → εx). Indeed, the crawling motion is very slow and the time variable needs to be "accelerated". Thus, we arrive at the following system of parabolic PDEs for a scalar phase field function ρ ε and the orientation vector P ε (describing the average orientation of the filaments that polymerize and push against the cell membrane): where λ ε (t) = 1 |Ω | Ω 1 ε 2 W (ρ ε ) + P ε · ∇ρ ε dx.
On the boundary ∂ Ω we impose the Neumann and the Dirichlet boundary conditions respectively ∂ ν ρ ε = 0 and P ε = 0, where Ω ⊂ R 2 is a bounded smooth domain.
Heuristically the system (1-3) describes the motion of cell membrane caused by the competition between mean curvature motion (due to surface tension) and the push of the cytoskeleton filaments against the membrane. The main issue is to determine the influence of this competition on the qualitative behavior of the sharp interface solution. The parameter β models this competition (as well as the adhesion to the substrate) which is why it plays a key role in the analysis of the system (1)(2)(3).
Equation (1) is a perturbation of the following Allen-Cahn equation [1]: This equation is a scalar version of the celebrated Ginzburg-Landau equation and it plays a fundamental role in mathematical modeling of phase transitions. It consists of the standard linear parabolic equation and a nonlinear lower order term, which is the derivative of a smooth double equal well potential Equation (4) was introduced to model the motion of the phase-antiphase boundary (interface) between two grains in a solid material. Analysis of (4) as ε → 0 leads to the asymptotic solution that takes values 0 and 1 in the domains corresponding to two phases separated by an interface of width of order ε, the so-called sharp interface. Furthermore, it was shown that this sharp interface exhibits mean curvature motion.
Recall that in this motion the normal component of the velocity of each point of the surface is equal to the mean curvature of the surface at this point. This motion has been extensively studied in the geometrical community (e.g., [18,20,17,4] and references therein). It also received significant attention in PDE literature. Specifically [9] and [11] established existence of global viscosity solutions (weak solutions) for the mean curvature flow. Mean curvature motion of the interface in the limit ε → 0 was formally derived in [21], [14] and then justified in [12] by using the viscosity solutions techniques. The limit ε → 0 was also studied for a stochastically perturbed Allen-Cahn equation (4) (e.g., in [22,27]). Note that equation (4) is closely related to another well-known model of phase separation, the so-called Cahn-Hilliard equation [6], which is a forth order reaction diffusion equation that models how two components of a binary fluid spontaneously separate and form domains of two pure fluids.
There are two distinguishing features in the problem (1)-(3): coupling and a nonlocal mass conservation constraint. We first comment on the coupling. The FitzHugh-Nagumo model, which is another prominent biological model, has a similar coupling feature, but in (1)-(3) the coupling occurs via spatial gradients (gradient coupling) of the unknown functions whereas in FitzHugh-Nagumo [30], [33] the two equations are coupled via lower order terms (unknown functions rather than their derivatives). There are several phase field models for the Allen-Cahn (also Cahn-Hilliard) equation coupled with another parabolic equation via lower order terms [12,7]. Our analysis shows that the gradient coupling results in novel mathematical features such as the following nonlinear nonlocal equation for the velocity of the interface curve Γ (t) derived below: Here V stands for the normal velocity of Γ (t) with respect to the inward normal, and κ denotes the curvature of Γ (t), c 0 is a constant determined by the potential W (c 0 = 3/2 for the specific choice (5) of the double-well potential), |Γ (t)| is the curve length, and the function Φ(V ) is given by (31). Next note that the term λ ε (t) given by (3) is a Lagrange multiplier responsible for the volume constraint (conservation of mass in the original physical problem [13]) Solutions of the stationary Allen-Cahn equation with the volume constraint were studied in [24] by Γ -convergence techniques applied to the stationary variational problem corresponding to (4). It was established that the Γ -limiting functional is the interface perimeter (curve length in 2D or surface area in higher dimensions). Subsequently in the work [31] an evolutionary reaction-diffusion equation with double-well potential and nonlocal term that describes the volume constraint was studied. The following asymptotic formula for evolution of the interface Γ in the form of volume preserving mean curvature flow was formally derived in [31]: Formula (7) was rigorously justified in the radially symmetric case in [5] and in the general case in [8].
Three main approaches to the study of asymptotic behavior (sharp interface limit) of solutions of phase field equations and systems have been developed.
When a comparison principle for solutions applies, a PDE approach based on viscosity solutions techniques was successfully used in [12,3,34,23,2,16]. This approach can not be applied to the system (1)-(3), because of the gradient coupling and nonlocal multiplier λ ε (t).
It is an open issue to introduce weak (e.g., viscosity or Brakke type) solutions in problems with constraints. Another technique used in such problems is Γ -convergence (see [32,22] and references therein). This also does not work for the system (1)-(3). The standard Allen-Cahn equation (4) is a gradient flow (in L 2 metric) with GL energy functional, which is why one can use the Γ -convergence approach. However, there is no energy functional such that problem (1)-(3) can be written as a gradient flow.
When none of the above elegant tools apply, one can use the technically hard approach of direct construction of an asymptotic expansion followed by its justification via energy bounds [26]. In Allen-Cahn type problems it typically requires a number of terms (e.g., at least five in [8]) in the expansion. Also only a special class of the initial data can be considered that contains all of the terms ("very well-prepared").
In this work we construct an asymptotic formula with only two terms: the leading one and the corrector (see e.g., (41)). The main ingredient is an appropriate choice of these terms which allows for good energy bounds for the corrector. Furthermore this choice leads to a reduction of the coupled system to a singularly perturbed nonlinear PDE (100) which for ε → 0 provides the sharp interface limit. This approach is rigorously justified in the 1D model problem, however we believe that this justification can be carried out in the 2D problem (1)- (3). For small β it is implemented via the contraction mapping principle; for large β it requires much more subtle stability analysis based on the spectral properties of a non self-adjoint non-local operator. Numerical simulations in the large β case shows non smooth velocity fields with hysteresis loops.
The main objectives of this work are: First, to show that there is no finite time blow up and that the sharp interface property of the initial data propagates in time. Theorem 1 establishes the existence of the solution to problem (1)-(3) on the time-interval [0, T ] for any T > 0 and ε sufficiently small ε < ε 0 (T ). It also shows that a sharp (width ε) interface at t = 0 remains sharp for t ∈ (0, T ). This is proved by combining a maximum principle with energy type bounds.
Second, to investigate how the presence of the gradient coupling and nonlocal volume constraints in (1-3) affect the sharp interface limit. We use formal asymptotics as ε → 0 to derive equation (6), which is the most important result for applications as it describes motion of the cell membrane.
It turned out that both the analysis and the behavior of the interface essentially depend on the parameter β which can be thought of as the strength of the gradient coupling. If β is small, then (1) and (6) can be viewed as a perturbation of Allen-Cahn equation (4) and curvature motion (7), respectively. In this case a suitable modification of techniques from [8] can be used to justify formal derivation of (6) and establish existence and uniqueness of its smooth solution.
On the other hand, if β is not small, then the standard techniques do not apply. Moreover, for large β equation (6) has multiple solutions because of non-monotonicity of the function c 0 V − β Φ(V ). Due to extreme technical complexity of the problem, a rigorous approach is developed for a 1D model problem. The derivation of the sharp interface equation is justified by developing special representations of the solution leading to an asymptotic reduction of the system to a single non-linear singularly perturbed PDE. In order to pass to the limit ε → 0 in this PDE, it is linearized around its steady states. Here the main technical tool is the spectral analysis of the linearized operator, which is non-local and non-self-adjoint. In particular, emergence of multiple solutions that may become discontinuous was established. Finally, numerical analysis of the original PDE system and of the sharp interface equation predicts the hysteresis phenomenon (which is not present in curvature motion with or without volume constraints because of linearity of these problems).
The paper is organized as follows. Section 2.1 is devoted to the well-posedness of the problem (1)- (2). In Section 2.2 the equation for the interface motion (6) is formally derived. Section 3 contains the rigorous justification of the sharp interface Limit in the context of a model one-dimensional problem.
2 Phase-field model of motion of cell. Well-posedness and formal derivation of sharp interface limit 2.1 Existence of the solution of (1)-(2) with ε-transition layer and no finite time blow up In this section we consider the boundary value problem (1)- (2). Introduce the following energy-type functionals The system (1)-(2) is supplied with initial data that satisfy: and The first condition (9) is a weakened form of a standard condition 0 ≤ ρ ε (x, 0) ≤ 1 for the phase field variable. If λ ε ≡ 0, then the maximum principle implies 0 ≤ ρ ε (x,t) ≤ 1 for t > 0. The presence of nontrivial λ ε leads to an "extended interval" for ρ ε 1 . The second condition (10) means that at t = 0 the function ρ ε has the structure of an "ε-transition layer", that is, the domain Ω consists of three subdomains: one where ρ ε ∼ 1 (inside the cell), another where ρ ε ∼ 0 (outside the cell), and they are separated by a transition layer of width ε (a diffusive interface). Furthermore, it can be shown that the magnitude of the orientation field P ε is small everywhere except the ε-transition layer (see (19)).
Remark 1 This theorem implies that if the initial data are well-prepared in the sense of (9)-(10), then for 0 < t < T the solution exists and has the structure of an εtransition layer. Moreover, the bound on initial data (9) remains true for t > 0. While it relies on a maximum principle argument, it also requires additional estimates on λ ε as seen from (13) below.
Proof: STEP 1. First multiply (1) by ∂ t ρ ε and integrate over Ω : Here we used the fact that due to (3) the integral of ∂ t ρ ε over Ω is zero and thus λ ε (t)∂ t ρ ε dx = 0.
Next, using the maximum principle in (1) we get (see Appendix A.3 for details): Let T ε > 0 be the maximal time such that and from now on we assume that t ≤ T ε .
Using (12), (14) and integration by parts we obtain We proceed by deriving an upper bound for the integral in the right hand side of (15). By (1) we have The following bounds are obtained by routine applying the Cauchy-Schwarz and Young's inequalities. For the the sum of the first three terms in (16) we get, Finally, in order to bound I 5 we first derive, then Thus, and using this inequality, (16) and (14) in (15), then substituting the resulting bound in (12) we obtain, for sufficiently small ε, STEP 2. Now we obtain a bound for the last two terms in (18). Taking the scalar product of (2) with 2kP ε + 4|P ε | 2 P ε , k > 0, integrating over Ω and using (14) we get We chose k := C 1 + 1 to obtain (18) and (19) we have the differential inequality, with a constant C > 0 independent of ε. Considering the bounds on the initial data and assuming that ε is sufficiently small, one can easily construct a bounded super-solutionG of (20) on [0, T ] such thatG(0) ≥ G ε . We now have, G ε ≤ C on [0, T ε ] for sufficiently small ε. By (13) and (17) we then conclude that T ε in (14) coincides actually with T when ε is small. The theorem is proved.

Remark 2
We briefly summarize the key points of the proof. Equation (1) is the Allen-Cahn PDE perturbed by a nonlinear quadratic term P ε · ∇ρ ε and the term λ ε which is constant in x. Therefore, the solution may exhibit a finite time blow up depending on the competition between Allen-Cahn part and the quadratic nonlinearity. The Allen-Cahn equation is a gradient flow corresponding to scalar Ginzburg-Landau energy. This energy admits well known estimates that rule out the blow up (by nonnegativity of energy). On the other hand even ODEs with quadratic nonlinearity are known to exhibit a finite time blow up. In presence of the second term blow up is ruled out by establishing differential inequality (20) (cf. standard Gronwall). The key difficulty here is due to the "fourth power" term 1 2 Ω |P ε | 2 |∇ρ ε | 2 dx in (12). Here a straightforward way would be to get L ∞ estimates on P ε directly, which is, however, not possible. Thus we first lower the power using integration by parts and maximum principle (14) to obtain a "third power term" in the RHS of (15). A subtle point here is that the maximum principle gives (13) which contains λ ε and one needs estimate (17) on λ ε to get bound (14) on an extended interval. Combining (15) with (1)- (2) leads to (20), which in turn rules out the blow up. (6) In this section we formally derive equation (6) for the 2D system (1)- (2). While the derivation is analogous to the single Allen-Cahn equation (e.g., [21], [8]), the gradient coupling in (1)-(2) results in a nonlinear term that modifies the mean curvature motion.

Formal derivation of the Sharp Interface Equation
Assume that that initial data ρ ε (x, 0) converge to the characteristic function of a smooth subdomain ω 0 ⊂ Ω . Next we want to describe the evolution of the interface where ω t is the support of lim ε→0 ρ ε (x,t). We will assume that the initial data coincide with initial values of asymptotic expansions for ρ ε and P ε to be constructed below.
Let X 0 (s,t) be a parametrization of Γ (t). In a vicinity of Γ (t) the parameters s and the signed distance r to Γ (t) will be used as local coordinates, so that where ν is an inward normal to Γ (t).
The inverse mapping to x = X(r, s,t) is given by is the limiting location of interface as ε → 0. Next we seek ρ ε and P ε in the following forms in local coordinates (r, s): Introduce asymptotic expansions in local coordinates: Now, substitute (22)-(24) into (1) and (2). Collecting terms with likewise powers of ε (ε −2 and ε −1 ) and equating them to zero we successively get: and the system (1)-(2) reduces to where κ(s,t) is the curvature of Γ 0 (t) and V 0 (t) := ∂ t r is the limiting velocity. The curvature κ appears in the equation when one rewrites the Laplace operator in (1) in local coordinates (r, s).
It is well-known that for every potential W (ρ) with two equal wells at ρ = a and ρ = b, a < b (in particular, for W (ρ) = 1/4ρ 2 (1 − ρ) 2 , a = 0, b = 1), there exists a standing wave solution θ 0 (z) of (25) which tends to b and a as z → ∞ and z → −∞, respectively. Moreover, all derivatives of the function θ 0 (z) exponentially decay to 0 as |z| → ∞ and θ 0 (z) is an eigenfunction of the linearized Allen-Cahn operator L u := −u + W (θ 0 )u corresponding to the eigenvalue 0. Then multiplication of (26) by θ 0 (z) and integration over z lead to the solvability condition for the right hand side of (26): Next we obtain the formula for λ 0 (t). It follows from (3) that Ω ∂ t ρ ε = 0. Substitute expansions (22) for ρ ε into Ω ∂ t ρ ε = 0 and take into account the fact that Thus, in order to satisfy the condition Ω ∂ t ρ ε = 0 to the leading order, V 0 (s,t) must satisfy V (s,t)| ∂ ∂ s X 0 (s,t)|ds = 0. Using this fact and integrating (28) with respect to s with the weight | ∂ ∂ s X 0 (s,t)|, we get Finally, the unique solution of (27) is given by The representation Ψ 0 (z, s,t) = ψ(z;V 0 (s,t))ν(s,t) yields where we have also taken into account the linearity of (30) in β . Now substitute (31) and (29) into equation (28) to conclude the derivation of sharp interface equation (6).

Rigorous asymptotic analysis for 1D model problem
The Sharp Interface Equation (6) formally derived in Subsection 2.2 exhibits striking non-uniqueness phenomenon for the normal velocity V . Namely, for large β the equation may have multiple roots. Note that combining together the curvature and integral (constant) terms in (6) yields the equation of the form (32) with Therefore, the uniqueness of the velocity field V that holds for small β is no longer guaranteed for arbitrary β . In this Section we analyze a 1D analogue of the original model which allows for rigorous analysis of the non-uniqueness phenomenon. The system is written as follows where F(t) ∈ C[0, T ] is a given function. We consider only bounded solutions: Analysis of the 1D problem (33)-(34) is a necessary step for understanding the original problem (1)- (2). Observe that motion of the interface in the 2D system (1)-(2) occurs in the normal direction, and therefore it is essentially one-dimensional. Thus, the 1D model (33)-(34) captures main features of the normal velocity in (1)- (2). However, the 1D model cannot contain the mass conservation constraint, because such a constraint in 1D is too restrictive and results in no motion. Therefore, the effects of curvature and mass conservation in (6) are modeled by a given function F(t). We believe that qualitative conclusions obtained for the 1D problem (33)- (34) apply for the 2D model (1)- (2).
with |v ε | 2 < c with c independent of ε and t.
The representation (37) has the leading term It has a transition region of "size slightly larger than ε", where it grows from 0 to 1 (see Fig. 1). The ε-interface x ε (t) can be viewed as a location ("center") of the transition region: time t. Remark 3 explains that a choice of ε-interface is not unique. For example, for a given ρ ε satisfying (37), every Our main goal is to prove the convergence as ε → 0 of the ε-interface x ε (t) to the sharp interface x 0 (t) whose velocity V 0 (t) =ẋ 0 (t) solves the sharp interface equation where Φ(V ) is a known nonlinear function given by (31). This equation can be formally obtained in the limit ε → 0 as in the Section 2.2.
Next we summarize key steps of our asymptotic analysis: 1. Choice of a special representation. The unknown function ρ ε is represented in the form (41), where θ 0 and ψ are known functions, and u ε is a new unknown function. Existence of the ε-interface x ε (t) and estimates on u ε uniform in ε and t are established in Section 3.2.
2. Reduction of the system to a single equation. The unknown function u ε is eliminated by proving that the third term in representation (41) is small. Next, we split P ε into two parts, P ε = A ε + B ε , where B ε depends on u ε but is small, and A ε does not depend on u ε . Thus, the original system (33)-(34) is reduced to Taking the limit ε → 0 in the system (39)-(40) is non-trivial because of the product term V ε (t)A ε .
3. Analysis of reduced problem. For sufficiently small β we prove that x ε (t) → x 0 (t) as ε → 0 by the contraction mapping principle. For larger β , system (39)-(40) further reduces to a singularly perturbed non-linear non-local equation. The limiting transition in this equation is based on the stability analysis of the semigroup generated by the linearized operator. This analysis leads to surprising results of non-uniqueness of velocities and hysteresis.

Asymptotic representation for ρ ε
In order to pass to the limit ε → 0 in (33)-(34) we further specify v ε in (37). Namely, we introduce the representation with new unknown function u ε satisfying and ψ ε (y,t) defined by , where ψ + and ψ − are solutions of the following ODEs with zero initial data. The purpose of the decomposition of the lower order term in (37) into two parts is two-fold.
Observe that the most interesting behavior of ρ ε occurs in the vicinity of the interface. So, roughly speaking, ideally we would like to "localize" our analysis by considering functions that are negligibly small outside the interface. However, the right hand side F(t) does not depend on spatial variable. Therefore, it does not decay at infinity which prevents ρ ε from being localized. This can be illustrated as follows: Thus, introduction of ψ ε leads to a "localized" behavior of new unknown function u ε which decays at infinity and, therefore, it allows one to work in Sobolev spaces on R.
In order to obtain bounds on u ε , the following Poincaré-type inequalities proven in Appendix A.1 are used: Here δ > 0 is arbitrary, k = 2, 3, 4, and w ε := u ε /(θ 0 ). Specifically, estimates for u ε are based on the identity: where L v := −v + W (θ 0 )v is the linearized Allen-Cahn operator that appears in the equation for u ε (see (51)).
Since θ 0 is a given function, we obtain Remark 3 The orthogonality condition (42) is achieved by an appropriate choice of ε-interface which can be explained as follows. Consider an ε-interface x ε (t) and substitute an ε 2 -perturbationx ε (t) = x ε (t) + ε 2 χ(t) into (41). Then re-arranging terms of order ε and choosing Thus, an ε 2 -perturbation of an ε-interface is an ε-interface in the sense of Definition 1. Moreover, for a given ε-interface one can find its ε 2 -perturbation that satisfies (42).

Reduction of the system to a single equation
The following theorem justifies the expansions (41) and will be used to obtain a reduced system for unknowns x ε (t) and P ε (y,t) by eliminating u ε .
Theorem 2 Let ρ ε , P ε be a solution of Problem (33)- (34) with initial data for ρ ε and P ε satisfying (36) and P ε ( and Then there exists an ε-interface Proof: , differentiate ρ ε from (41) in time, multiply by θ 0 (y) and integrate in y over R. Then we get (u ε is absent thanks to (42)) Sinceẋ (50) is an ODE for the unknown function x ε (t). Standard theorems from ODE theory imply existence and uniqueness of x ε (t). STEP 2. (equation for u ε and energy estimates for terms with u ε ) Substitute expansion (41) into equation (33): Introduce the function w ε satisfying u ε = θ 0 w ε . Multiply equation (51) by u ε and integrate in y over R.
Here we used Poincaré inequalities (see (160) and (162) in Appendix A.1) to bound terms with |u ε | 3 and u 4 ε : To estimate the term with ψ we use boundedness of ψ where δ is a small number which does not depend on ε. Thus, we have We will estimate terms containing P ε in the next step. STEP 3. (representation for P ε and energy bounds for terms with P ε ) To estimate terms with P ε we introduce the representation equipped with initial conditions A ε (y, 0) = P ε (y, 0) and B ε (y, 0) = 0. The maximum principle implies that |A ε (y,t)| ≤ max y∈R |p ε (y)| + |β | max y∈R |θ 0 (y)|. Multiplying equation (56) by B ε and integrating in y over R we get Then terms with P ε are estimated as follows: By adding (57) to (54) we get: STEP 4. (boundness of u ε ) Inequality (58) implies where V (t) := u 2 ε dy. By taking t such that V (t) ≤ 2C 1 + V (0) for t ∈ [0,t ) and εC 3 + ε 3 C 4 < C 1 /2 we get: Thus, for a given T we can take t = T and obtain the desired L 2 -bound, uniform in t, on u ε : The theorem is proved. We now derive the equation which in the limit ε → 0 becomes the desired sharp interface equation (38). Multiply (51) by θ 0 (y) and integrate in y over R to obtain whereÕ Terms with ∂ t u ε and u ε −W (θ 0 )u ε vanish due to the following identities: Proposition 1 Let ρ ε and P ε be solutions of (33-36). Then O ε (t) introduced in (60) satisfies the following estimate This follows from pointwise bound ρ ε L ∞ < C (obtained by maximum principle applied to (33)) and representation u ε = ε −1 (ρ ε − θ 0 ) − ψ ε . Using inequality (61) and Theorem 2 we obtain a bound for the first integral term in O ε (t): This implies Integrating in t we get STEP 3. It remains to bound the following terms in O ε (t): A ε ψ ε θ 0 dy + A ε u ε θ 0 dy + B ε ψ ε θ 0 dy + B ε u ε θ 0 dy + B ε (θ 0 ) 2 dy. (63) Applying (46) we bound from above the first and the third term by a constant C. Integrating the second term by parts and using the maximum principle for A ε (max Integration by parts of the fourth term and taking into account (62) yields From (57) we derive a uniform bound for B ε 2 L 2 (R) < C which in turns implies that Thus, all five terms in (63) are bounded and the proposition is proved.
Thus, we reduced the original system (33)- (34) to whereÕ ε (t) and O ε (t) are bounded in L ∞ (0, T ). Observe that taking the formal limit ε → 0 in the system (64,65) leads to a system of two equations equivalent to the sharp interface equation (38). Indeed, the formal limit ε → 0 in (65) is nothing but (30) whose unique solution is ψ(y;V (t)). Then substituting this function into the limit of (64) yields (38).

Sharp
Interface Limit for small β by contraction mapping principle.
In the following Theorem that establishes the sharp interface limit for sufficiently small β we assume that initial data P ε (εy, 0) = A ε (y, 0) are bounded in L 2 by a constant C independent of ε: Theorem 3 Let (A ε ,V ε ) be solution of the reduced system (64)-(65) withÕ ε , O ε ∈ L ∞ (0, T ), and initial conditions satisfying (66). Then there exists β 0 > 0 such that for where V 0 is the unique solution of (38), and β 0 depends only on θ 0 (y) which is given explicitly.
Proof: STEP 1 (construction of boundary layer at t = 0). Introduce the function η ε (y,t) that solves the following system (it will be shown below that η ε is the boundary layer) Multiplying equation (68) by η ε and integrating in y we get Integrating the equality (70) in t we obtain STEP 2 (Resetting of (64)-(65) as a fixed point problem). Define F ε : where A is the unique solution of Here V bl := 1 c 0 +εÕ ε (θ 0 ) 2 η ε dy and it accounts for the boundary layer at t = 0, and by (71) It follows from the construction of F ε that V ε − V bl is a fixed point of this mapping. Next we prove that, for sufficiently small β , F ε is a contraction mapping. Consider V 1 ,V 2 ∈ L ∞ (0, T ) and let A 1 , A 2 be solutions of (73)-(74) with V = V 1 and V = V 2 , respectively. The functionĀ := A 1 − A 2 solves the following problem Multiplying equation (76) byĀ and integrating in y we get In order to estimate A 2 2 dy note that for every solution A of (73) the following bound holds Indeed, multiplying (73) by A and integrating in y we get: This equality implies Integrating (81) and taking into account that |A(y, 0)| 2 dy = |ψ(y;V 0 (0))| 2 dy < Cβ 2 we obtain (79). Using (79) in (78) we get which, due to zero initial conditions forĀ, implies Finally, to prove that F ε is a contraction mapping we use (83): Choosing β sufficiently small, q := √ Cβ < 1 we get that F ε is a contraction mapping. As follows from the derivation of (79), the constant C depends only on θ 0 (y).

STEP 3. Observe that
Thus, Thus, it is left to prove that This will be done in the next step.
(88) Since the approximation V ε 0 no longer solves the sharp interface equation (38), we introduce a function F ε (t) by Then V ε 0 solves (38) with F ε (t) in place of F(t). The following convergence is straightforward Let A be the solution of (73)-(74) with V = V 0 (t). Then D ε (y,t) := A(y,t) − ψ(y;V ε 0 (t)) satisfies the following equality Integrating this differential inequality we get Taking into account the facts that Finally, we establish the L ∞ -convergence of F ε (V 0 ) to V 0 : (95) The convergence of the first term follows from (90) and (94). The convergence of the second term directly follows from (87).
Thus, (86) is established and the theorem is proved.
Remark 4 Theorem 3 shows uniform convergence of V ε to V 0 on [δ , T ] for δ > 0 and convergence A ε (y,t) to A 0 (y,t) := ψ(y;V 0 (t)) in C([δ , T ]; L 2 (R)). This result cannot be improved to uniform convergence on the whole interval [0, T ] because the only condition on the initial data A ε (y, 0) is (66), and A ε (y, 0) does not necessarily converge to A 0 (y, 0) in any sense. This shows that the solution A ε (y,t) for small t exhibits a boundary layer behavior. This behavior is captured by η ε given by (68)-(69).

Reduction to a stability problem
For larger β the contraction principle no longer applies and both analysis and the results become much more complex. In particular, surprising features of the motion of the interface such as discontinuities of the interface velocity and hysteresis will be established. Here the stability analysis of the semigroup generated by a non-local non self-adjoint operator is used in place of the contraction mapping principle.
In the case where β is not small, solutions of (38) are no longer unique. However, the original PDE problem (33)-(34) (as well as the reduced system (64)-(65)) has the unique solution. This indicates that analysis for large β must be complemented by a criterium of how to select the limiting solution of equation (38) among all solutions of this equation.
As a first step, we neglect termsÕ ε (t) and O ε (t) in the reduced system (64)-(65) and study the following system: (97) In the system (96)-(97), f ε replaces A ε from (64),(65). Note that a formal limit ε → 0 in the system (96)-(97) gives Substituting (96) into (97), we rewrite the system (96)-(97) as a single equation: In the limit ε → 0 this equation reduces to single PDE with unknown f 0 : Taking the formal limit is justified below for passing from (100) to (101). Equation (100) is a singular perturbation of (101) and both equations are nonautonomous. It is well-known that singular limit problems, including non-autonomous equations, can be reduced to the analysis of large time behavior of autonomous equations. To illustrate this, recall a standard example of an ODE with a small parameter ε from [19]: Assume that there exists the unique root φ (t) of F , i.e., 0 = F (φ (t),t), t ∈ [0, T ].
Then the singular limit φ (t) = lim ε→0 z ε (t) holds provided that φ (t) is a stable root, i.e., all solutions u(τ) of an autonomous problem du(τ) dτ = F (u(τ),t) (t is fixed) converge to the large time limit φ (t): lim τ→∞ u(τ) = φ (t). Note that the problem (102) has two time scales: a slow time t and a fast time τ. Also the large-time limit corresponds to τ → ∞ for a fixed parameter t.
Note that the equivalence of singular and large-time limits is straightforward for the singularly perturbed autonomous problems (F does not depend on t in (102)). In this case, the simple rescaling τ := t/ε u(τ) := z ε (ετ) reduces the singular limit problem to a problem of stability of steady state. We use a similar approach to justify the transition as ε → 0 (singular limit) from non-autonomous equation (100) to equation (101) by introducing three time scales: slow, fast, and intermediate. More precisely, we employ the following three step procedure: (i) partition the interval [0, T ] by segments of length √ ε on which the equation (100) is "almost" autonomous (F(t) is "almost" constant on each of these intervals); (ii) on the first interval (0, √ ε), by appropriate scaling τ = t/ε and stability analysis find large-time asymptotics τ → ∞ (here we used equivalence of singular and large-time limits for autonomous equations); (iii) use the asymptotics found in (ii) as initial conditions for the next interval ( √ ε, 2 √ ε), repeat step (ii) on this interval, and continue to obtain global asymptotics on [0, T ]. A crucial ingredient here is an exponential stability of the linearized problem which prevents accumulation of error due to the bound (126) from Lemma 1.

Stability analysis of the rescaled equation
In this subsection we focus on stability as τ → ∞ of the following equation obtained as described in step (ii) above by rescaling fast time τ = t/ε in the unknown f ε in (100) and "freezing" time t in F(t): Here t ∈ [0, T ] is considered to be a fixed parameter. We describe the large-time behavior (i.e., as τ → ∞) of solutions of (103) using both numerical simulations and rigorous analysis. In particular, we show that solutions of (103) converge to a steady state that solves (101). As follows from the discussion in the previous subsection, this convergence is necessary to pass to the singular limit ε → 0 in the equation (100). It will be shown below that the equation (101) has a finite number of solutions which varies depending on parameters β and t. We will provide a criterion which determines whether a given steady state of the problem (103) is stable or unstable. The full description of large-time behavior of (103) is first obtained numerically, and then it is justified by linear stability analysis for certain values of V 0 (V 0 is defined in (98)).
For a given solution f of (103) we define its velocity as follows This equality coincides with the equation (96) for velocity V ε (t).
The steady state velocity V st corresponding to a steady state f st is determined by Using (31) and (101), this definition of V st can be rewritten as follows This equation is nothing but the sharp interface equation (38). Once V st is found from (105), the steady state f st is uniquely determined by the equality f st (y) = ψ(y;V st ). One can show that Φ(V ) is a real analytic function of V decaying (with all its derivatives) as |V | → ∞. It follows then that (105) has at least one and at most finite number of solutions.
We begin with two numerical observations. First, depending on parameters β and t, equation (105) possesses one, two or three roots. Numerical simulations show that Φ(V ) is an even function monotonically decreasing to 0 as |V | → ∞, see Fig. 2 (left). For small β the function c 0 V − β Φ(V ) is monotone and invertible as a function of V . Thus, in this case, (105) has only one root. For larger β , function c 0 V −β Φ(V ) has two local extrema, local maximum V max and local minimum V min . In this case, Here The case of larger β , when the number of roots varies depending on F(t), is illustrated by Fig. 2 (right). Second, for every solution f of (103), and for arbitrary β and F(t), the velocity V (τ) defined by (104) converges to a solution V st of (105) as τ → ∞. Moreover, V st satisfies the following inequality: To illustrate this observation consider the case when (105) has three distinct roots V 1 < V 2 < V 3 (see Fig. 3, left). In numerical simulations we observe that steady state velocities V 1 and V 3 are stable, and V 2 is unstable. Specifically, a velocity V (τ) given by (104) converges to either V 1 or V 3 for various initial conditions (see Fig. 3, right). Note that V 1 and V 3 satisfy (107).
Note that our primary concern is not the convergence of velocity V (τ) but rather the convergence of f (y, τ) that solves equation (103) to a steady state.
Next we study stability of steady states by means of spectral analysis of the linearized operator. Let f st be a solution of (101) (a steady state of (103)) and V st be the corresponding steady state velocity defined by (105). The linearization of equation (103) around f st has the following form: where T (V ) : L 2 (R) → L 2 (R) is a linear operator parameterized by V ∈ R and defined by where ψ (y;V ) is defined in (30). Operator T (V ) is a perturbation of a local operator A (V )g := −g − V g + g by a non-local rank one operator P( The spectrum σ (A (V )) of operator A (V ) is described by the following straightforward proposition.

Proposition 2
The spectrum σ (A (V )) consists only of its essential part: Next we describe the spectrum σ (T (V )) of operator T (V ).
Theorem 4 Consider the part σ p (T (V )) of the spectrum σ (T (V )) laying in C\σ ess (A ). Then σ p (T (V )) is described by Moreover, all λ from σ p (T (V )) are simple eigenvalues.

Proof:
We suppress dependence of A , T , and η on V for the sake of brevity.

Proof:
Again for the sake of convenience, we suppress dependence of A , T , and η on V . STEP 1. We prove the following equality which is used in STEP 2.
We next derive stability conditions for the linearized equation (108). Consider λ = λ r + iλ i . Denote by H λ (k) the integrand in (113) and rewrite it in the form , then all eigenvalues λ ∈ σ p (T (V )) have positive real part, Reλ > 0.
(ii) There existsV < 0 such that for all V <V all eigenvalues of T (V ) have positive real part.
Remark 6 Condition V > − √ 2 is a technical assumption in the proof which guarantees that integral (120) is negative. However, numerical simulations show that integral (120) is actually negative for all V .
which is negative or zero and, thus, cannot be equal to 1, so equality (110) does not hold and, in particular, there does not exist eigenvalues with negative real part. Thus, part (i) is proved.

Main result for 1D interface limit
In this subsection we formulate the main result on the 1D sharp interface limit in the system (33)-(34) for arbitrary β and illustrate the qualitative behavior of the limiting velocity V 0 (t) = lim ε→0 V ε (t).
To formulate the main result we need the following definition.
Definition 2 Define the set of stable velocities S D by where T (V ) is the linearized operator given by (109).

Remark 7
We summarize here what is known about the set S D of stable velocities. According to Proposition 3 the set S D is contained in Moreover, based on numerical computations of integral (111) we conjecture that Next theorem is the main result for 1D sharp interface limit. Introduce the following conditions: (C1) (Conditions on solutions of the sharp interface equation) Let V 0 ∈ S D be a solution of V 0 − β Φ(V 0 ) = F(0) and let [0, T ] be a time interval such that ∃V (t) ∈ S D a continuous solution of (C2) (Well-prepared initial conditions for P ε ): P ε (εy, 0) = p ε (y) and p ε −ψ(·;V 0 ) L 2 < δ with a small constant δ > 0 independent of ε. The function ψ = ψ(y;V 0 ) is defined by (30). Theorem 5 justifies the sharp interface equation (122) for ANY β . Its proof consists of two steps: (i) reduction to a single equation (nonlinear, singularly perturbed) which is done in Section 3.2 and (ii) passage to the limit in this equation based on stability analysis presented below, which is the main ingredient of the proof.
Proof: of Theorem 5. Rewrite (64)-(65) in the form of a single PDE where κ (i) ε , i = 1, 2, can be computed from (64)-(65) and satisfy the following bounds (actual structure of κ We next pass to the limit in equation (123) using exponential stability (established in (132)) of the semigroup corresponding to the linearized operator. The following local stability result plays the crucial role in the proof.
This Lemma shows that if initial data are at distance at most δ from ψ (in the L 2norm), then the solution A ε (y,t + ετ) approaches ψ exponentially fast in τ (first term in the RHS of (125)) with a deviation that is bounded from above independently of t (described by the second and the third terms in the RHS of (125)). The conclusion of Theorem 5 immediately follows from this Lemma. Indeed, consider the time interval (0,t 1 := √ ε). Then by Lemma 1 we obtain where m denotes the modulus of continuity of F on [0, T ]. Choose ε small enough so that the RHS of (126) admits the bounds Similarly for intervals ( To complete the proof of Theorem 5 we again use Lemma 1 to bound A ε (·,t) − ψ(·;V (t)) 2 L 2 for t ∈ (t i ,t i+1 ). Proof of Lemma 1. Consider the function f ε (y, τ) := A ε (y,t + ετ) − ψ(y;V (t)). It follows from (123) and (30) that f ε satisfies the following PDE where T = T (V (t)) andκ ε | ≤ C with a constant C independent of t and ε. Then, introducing the semigroup operator e −T τ , τ > 0 in L 2 (R), we have by Duhamel's principle: In order to proceed with the proof of Lemma 1 we need to prove the exponential stability of the semigroup e −T t and establish its consequences in the following Lemma 2 There exists ω > 0 such that the following inequality holds with constant C independent of g. Moreover, both constants ω and C can be chosen independently of t (recall that T = T (V (t)) depends on t).
In a standard way, multiplying (137) by v and integrating in y we get Then an application of Gronwall's inequality yields the uniform bound v(·, τ) L 2 ≤ C for 0 ≤ τ ≤ 1.
Using this bound in (139) we derive It follows that v (·, τ 0 ) ≤ C g for some 0 < τ 0 ≤ 1. Then by the semigroup property: where we have used (132). STEP 3. Now using (136) and the Cauchy-Schwarz inequality we get To complete the proof of (131) consider It follows from the definition off thatf (y) =ṽ(1, y), where v solves In the standard way, multiply equation (146) by v and integrate in y to obtain It follows that ṽ(·, 1) 2 and τ τ−1 e −T (τ−τ ) g (τ )dτ 2 L 2 = ṽ(·, 1) 2 Combining (149) with (144) completes the proof of the Lemma 2. Now we apply Lemma 2 to (130) to obtain the bound where δ ε = m 2 (εs) + ε 2 and m denotes the modulus of continuity of F on [0, T * ]. Let α(τ) be the right hand side of (150), theṅ We can choose small enough δ that if f ε (·, 0) ≤ δ , then for sufficiently small ε the functionᾱ(τ) being the solution of the ODĖ Note thatᾱ(τ) is given bȳ Recall that α(0) ≥ 0, δ ε ≥ 0 and δ ε (τ) → 0 as ε → 0. Therefore, if 0 ≤ α(0) ≤ Cδ , then for sufficiently small ε (such that δ ε ≤ δ ) we have 0 ≤ᾱ(τ) ≤ C δ for all τ ∈ [0, ε − 1 r ]. It follows that for sufficiently small δ α ≥ −ωᾱ +C(δ ε +ᾱ)ᾱ + δ ε , and consequently (151) does hold. To conclude the proof of the Lemma 1 (and, thus, to establish Theorem 5) it remains to note that In view of the above analysis the equation (38) for large β may have many solutions of quite complicated structure (e.g., discontinuous). Therefore, we need to introduce a criterion for selection of the "correct" solutions that are limiting solutions Fig. 4 (SKETCH) Left: Until F < F max we stay on the same the left branch, When F exceeds F max we jump on the right branch; Right: Until F > F min we stay on the same the right branch; When F exceeds F min we jump on the left branch. 0 ≤ t 0 < t 1 < t 2 < t 3 ≤ T ; F max and F min are defined in (106) to the problem with ε > 0. This is analogous, e.g. to viscosity solutions of Allen-Cahn when physical solutions are obtained (by regularization) in the sharp interface limit ε → 0, [12].
We now introduce such a criterion based on numerical observations and suggested by the stability analysis depicted on Fig. 4. Define the left velocity interval B L := (−∞,V max ] and the right velocity interval B R := [V min , ∞) for stable velocities V .
Assume for simplicity of presentation that function F(t) ∈ C[0, T ] is strictly increasing. Then the solution of (38) is chosen based on the following two criteria (Cr1) if V (0) ∈ B R , there is a unique V (t) ∈ B R satisfying (38) for all t ∈ [0, T ]. Note that this V (t) is the only solution which is continuous and never enters the "forbidden" interval [V max ,V min ] (Cr2) if V (0) ∈ B L , then for any t ∈ [0, T ] the solution V (t) of (38) is chosen in the left velocity interval B L , unless it is impossible (F(t) > F max , where F max is defined in (106)). In the latter case V (t) is chosen from the right velocity interval B R .
Intuitively, evolution of the sharp interface velocity can be described as follows. Consider for example the left part of Figure 4. As time evolves, the velocity increases along the left green branch until it reaches V max , then it jumps (along the horizontal red dashed line) to the solution of (38) on the right green branch, and continues increasing along this branch.
Consider first y > 0, we have e −κy |v(y)| = |κ The same arguments for y < 0 complete the proof the Theorem. Treating similarly y < 0 and the using (154) one derives (162).