{'fraction_non_alphanumeric': 0.04724991436967064, 'fraction_numerical': 0.016783544554812426, 'mean_word_length': 4.164029044870602, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 8, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 25, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Revisiting Shikhmurzaev's approach to the contact line problem

In this paper, we revisit a model for the contact line problem which has been proposed by Shikhmurzaev (1993). In the first part, in addition to rederiving the model, we study in detail the assumptions required to obtain the isothermal limit of the model. We also derive in this paper several lubrication approximation models, based on Shikhmurzaev's approach. The first two lubrication models describe thin film flow of incompressible fluids on solid substrates, based on different orders of magnitude of the slip length parameter. The third lubrication model describes a meniscus formation where a wedge-shaped solid immerses in a thin film of fluid.


Introduction
In this paper, we revisit Shikhmurzaev's theory for fluid interfaces [25]. In particular, we discuss in detail the conditions required to be able to assume that the temperature of the system is constant. We derive also some thin film approximations of the model.
To explain the background briefly, a 'contact line' is formed when three distinct phases of matter meet. If the phases are all fluid (liquid or gas) phases, that is, two immiscible fluids are in contact with another liquid or a gas, this contact line is freely deformable in space, while it is bound to move on a given surface, if one of the phases is a solid. We are interested in this latter case. The process of displacement of one fluid by another immiscible fluid from a solid surface in the framework of classical fluid mechanics is known as the 'moving contact line' phenomenon. This contact line motion appears in many natural and technological processes, for example, coating of solids by liquid films, photographic films or polymer processing.
The mathematical modelling of the dynamic moving contact line encounters a fundamental difficulty, compared to the equilibrium state. The typical no-slip boundary condition for velocities of both the fluids on the solid surface gives rise to a non-integrable singularity of the shear stress, that is, to an infinite force exerted by the fluid on the solid (cf. [21], [12]) which is considered not to be physical. The pressure of the fluid at the contact line becomes singular as well. This problem is known in the literature as the "moving contact line paradox".
Being important from both theoretical and practical point of view, this problem has been investigated from different perspectives. In order to have a well-defined problem in the presence of a contact line, one must address two fundamental questions: (1) how to determine the (dynamic) contact angle i.e. the angle between the free surface and the solid boundary as the contact line moves, and (2) how to remove the divergence of the dissipation of the energy formula. In determining the contact angle, it is important to take into account that the contact angle is not merely a material property, but it depends also on the velocity of the contact line (as suggested by experiments). The survey of Dussan [11] discusses different approaches by physicists and chemists in determining the dynamics of a contact line. It is broadly accepted that the dynamic and equilibrium contact angles θ d and θ e are related via v c = F (cos(θ d ) − cos(θ e )), (1.1) where v c is the contact line velocity (along the solid), and F is some increasing function such that F (0) = 0. But there have been different suggestions for the precise form of F in equation (1.1). Blake-Haynes [6] derived (1.1) by using thermodynamic and molecular kinetics while Cox [9] used matched asymptotic analysis and hydrodynamic arguments. Also Ren-E [24] derived (1.1) from thermodynamic principles applied to constitutive equations. On the other hand, to remove the shear-stress singularity near the contact line, a widely used approach is to consider the Navier boundary condition which says that the slip of the fluid along the solid is proportional to the tangential stress acting from the fluid on the fluid-solid interface, i.e.
v L · n = v S · n, τ · S · n = β v L − v S · τ, with some parameter β > 0 which is known as the friction parameter, or the slip coefficient. Variants of the above slip condition concerning different assumptions on β have been considered in a number of works, for example in [19,20,17,1]. The above references are by no means exhaustive and we refer to [26,Section 3.4] for further reading. Another model has been introduced by Y. Shikhmurzaev in [25], based on the approach proposed by Bedeaux et al [2], to generalize the methods of irreversible thermodynamics to systems containing interfaces. The basic idea of Shikhmurzaev's theory is to consider a three phase model (cf. Figure 1) and describe the interfaces as surfaces of zero thickness, separating the bulk phases. Although the thickness of these interfacial layers is negligible compared to the characteristic length scale of the bulk flow, they can possess some intrinsic surface properties (e.g. surface tension) which play an important role in the overall dynamics via capillary effects. This approach formulates in a systematic way the boundary conditions at the interfaces for the equations describing the bulk phases. Similarly, the contact line can be defined as boundaries of the interfaces. Incorporating the fundamental laws of fluid mechanics and the theory of irreversible thermodynamics, in addition to the classical conservation equations of mass, momentum and energy, different constitutive equations and transport coefficients are obtained, giving the required boundary conditions. In order to model the thermodynamics of the interfaces, the idea of Onsager's principle has been used, as for the bulk phases. As a consequence, the (Navier) slip boundary condition is deduced at the fluid-solid interface which makes the force on the solid surface bounded, solving the paradoxical phenomena described in the beginning. Note that this derivation also makes clear that the slip boundary condition is not just an empirical assumption as often thought of, but rather can be formulated rigorously. Secondly, the dynamic contact angle is determined as part of the solution, given by a generalization of Young's law. The main novelty of this model is that it incorporates the aspect of interface formation and disappearance at the contact line, apart from addressing the two fundamental issues of contact line problem. This model is very appealing and has a great potential to explain several physical phenomena such as wetting, coalescence, cusp formation etc (cf. the introduction of [4]).
In Section 2, we rederive the model in situations which might contain contact lines at the boundary of the solid, liquid and gas phases. We made the derivation as mathematically rigorous as possible, using the physical ideas of Bedeaux et al [2], as well as Shikhmurzaev's. One of the issues that we try to clarify is the boundary conditions at the contact lines. It was indicated in [4] that there is a missing boundary condition which is not included in the original derivation of Shikhmurzaev (in [25], [26]). In [7], the authors also indicated the requirement of an additional boundary condition, although the one derived in [7] is slightly different from the one in [4]. In our derivation, we obtain the need for the additional boundary condition and we obtain the same one suggested in [7]. This additional condition is required in order to obtain a well-posed mathematical problem.
In Section 3, we derive different lubrication models using Shikhmurzaev's approach to the contact line problem. Here the motivation is to obtain in the particular setting of thin films, a simpler reduced model than the original one, in which the differences between this approach and the other solution models of the contact line problem (for example, Navier boundary condition, disjoining pressure, precursor film etc) can be compared. Thin film equations have been discussed in many works, both experimentally and theoretically, taking into account different theories and approaches, and to answer different questions. Derivation of the classical thin film equation in case of a droplet can be found in, for example [17,21], while in [15], the authors have considered thin film equations introducing the slip parameter in the model. On the other hand, there are other lubrication models, such as of a thin fluid in the Hele-Shaw cell [8], of Darcy flow [16], rupture of thin film subject to van der Waals forces [13], considering different orders of magnitude of the slip-length parameter [23], for an active gel [22], for surfactant driven viscous film [14]. In [27], the author has discussed a rupture phenomenon in the framework of the thin film approximation. Also an approximation of Shikhmurzaev's model for a gravity driven flow on an inclined plane has been considered in [5]. In particular, the asymptotic behaviour of the solutions in a thin film setting has been computed. We refer to the review paper [30] for an overview. The classical work showing the existence of a weak solution corresponding to the zero contact angle goes back to [3]. Also we mention the recent work [18] where the authors establish the existence of a global in time solution of the full contact line problem, provided the data is sufficiently close to the equilibrium value and the equilibrium contact angle may vary between 0 and π.

Setting
In the following, we summarize the main ideas in the derivation of Shikhmurzaev's model. Specifically, we reformulate the model in the notations that we will use in rest of the paper and discuss in detail some of the assumptions, in particular the isothermal limit, that has been used in most of Shikhmurzaev's work. Let us consider a simple model consisting of a liquid, gas and solid phase (more general situations with two liquids can be found in [26]). The model includes the evolution of interfaces separating two bulk phases.
Let us consider a liquid with density L . As we assume that our system can be described as a continuous medium, the properties of the fluid can be completely characterized by the velocity v L and the temperature T at each point (x, t). Concerning the gas, it has density G and pressure p G at each point, while for the solid body, we assume that it is union of a finite number of rigid pieces and therefore we do not need to describe its deformations. Nevertheless, one could consider situations where there is a relative motion between the different parts that composes the solid boundary, for instance a container whose part of the boundary is a mobile piston. The rigid solid has a temperature T , as well as density S .
The interaction between the different phases takes place through the liquid-gas, liquidsolid and gas-solid interfaces which are denoted by Σ 1 , Σ 2 , Σ 3 respectively. The three interfaces meet at a line C which is termed as the contact line (cf. Figure 1). We assume in the following that all the interfaces Σ i are contained in smooth surfaces and that the contact line C is a smooth curve. More precisely, we assume that the gas and the liquid fill a domain Ω with smooth boundary ∂Ω. Then ∂Ω = Σ 2 ∪ Σ 3 . Note that, in general, the domain Ω changes in time. This change of Ω happens if the domain translates rigidly in space or if Ω changes its shape, or a combination of translation and deformation of Ω. We can also have the situation where Ω does not change in time, but ∂Ω has a non-zero velocity. This happens for instance in the case in which the boundary of the container is a rotating cylinder. We characterize in the following the motion of the solid domain by a velocity function v S defined on ∂Ω. We will assume in the following, The parameter σ characterizes uniquely each material point of the boundary. Then the velocity of each boundary point is given by v S (σ, t) = ∂Γ ∂t (σ, t). Notice that the assumption that the boundary of the container consists of rigid parts imposes stringent constraints on Γ(σ, t) and v S (σ, t). We remark also that the domain Ω does not change in time if v S · n = 0 where n is the normal vector to the boundary ∂Ω.
These interfaces Σ i can be thought of as two-dimensional continuous media that can be described by a set of macroscopic variables. In our problem, the variables that we use to characterize these interfaces are the local surface densities s i , i = 1, 2, 3, the velocities v s i , i = 1, 2, 3 and the temperature T . We emphasize that the surface densities s i are not mass per unit volume, but mass per unit surface. Notice that we expect the mass for the interfaces to be small compared to the bulk, i.e. if denotes the thickness of the interface, then s i L . Moreover, denoting by L the characteristic length scale of the full system, the mass in the liquid is given by m L = L L 3 while the mass of the interface is m s 1. However, the interfaces are able to produce relevant mechanical forces (for example, surface tension) that can be strongly dependent on the densities s i , which play an important role in the properties and dynamics of the interfaces. This is the reason that we keep s i in the model. A relevant feature of Shikhmurzaev's approach is that s i can take different values at different points of the interfaces. Concerning the contact line, the velocity is denoted by v c , defined as in (A.4). At this point, we assume that the system is in local equilibrium, i.e. all the thermodynamic relations (such as Gibbs relation) hold for every infinitesimal region of the space both in the bulk and at the interfaces.

Assumptions
The main assumptions that we will use in the derivation are the following. We restrict ourselves to the case of an incompressible liquid, i.e. the liquid density L is assumed to be constant (the compressible case could be considered similarly).
Since the gas density is much lower than that of the liquid, we would like to consider G = 0. On the other hand, we want to have a gas with pressure comparable to that of the liquid. These gas properties have to be understood in a suitable asymptotic sense and the precise meaning and the consistency of these properties is explained in detail later in Section 2.2. We remark here that the assumption G L leads us to assume that both the internal and kinetic energy of the gas per unit of volume are negligible compared to the analogous quantities for the liquid. We will also see later that the heat transport in the gas is much smaller than in the liquid for similar differences of temperature.
The density of the solid body S is also assumed to be constant. As the solid is considered rigid, its velocity v S can be characterized by the translational velocity U(t) and rotational velocity ω(t), We also assume that the temperature across all the phases, the bulk, at the interface and at the contact line, is continuous, but not necessarily constant. That is, for We further assume that the total mass of the contact line is negligible, and also that it contains neither momentum nor energy. However, suitable conditions concerning the conservation of mass, momentum and energy hold on it. This is not the only possibility; In the re-formulation of Shikhmurzaev's model in [7], the contact line density has been assumed to be non-zero.
In the Appendix, we precise the notations used in this work and recall the necessary definitions and formulas from differential geometry. Throughout this article, we use the superscripts L, S, G, s to denote the quantities in the liquid, solid, gas and at the interfaces respectively.

Conservation laws
With the previous assumptions, the total mass density of the system becomes, where χ L , χ S are the characteristic functions of the region filled by the liquid and the solid respectively and δ Σ i denotes Dirac measures, supported on the surfaces Σ i , i = 1, 2, 3. For later use, we also introduce the characteristic function χ G of the region occupied by the gas. Notice that ∂(supp χ L ) = Σ 1 ∪ Σ 2 . Similarly, the total momentum density of the system is given by, Note that in quantities (2.3) and (2.4), the contribution of the gas to the mass and momentum density has been completely neglected, due to the above assumptions in Section 2.1.1. Further, the total energy density is Here e ξ , e s i are the total energies per unit mass, 1 2 |v ξ | 2 , 1 2 |v s i | 2 are the kinetic energies per unit mass and u ξ , u s i are the internal energies per unit mass, for the bulk and the interfaces respectively.
We neglect in (2.5) the internal energy of the gas compared with the one of the liquid. Indeed, for the changes of internal energy per unit mass, we have where c G m , c L m are the specific heat capacities of the gas and the liquid respectively. To obtain the above quantities for unit volume, we must multiply with the densities of the gas and the liquid respectively. Since the specific heat capacities of gas and liquid are of the same order of magnitude 1 , it follows that G u G L u L . Also, the kinetic energy per unit volume of the gas is much smaller compared to the liquid if v G , v L are comparable, i.e. 1 2 G |v G | 2 1 2 L |v L | 2 . Therefore, we get that G e G L e L and thus, we ignore the contribution of the gas to the energy and hence to the entropy.
We now derive the evolution equations for our system. These are obtained from the following basic conservation equations for mass, momentum and energy: Here, the momentum flux is defined as, 10) and the total energy flux as, The terms ξ v ξ ⊗ v ξ , s i v s i ⊗ v s i denote the transport of momentum in the bulk and at the interfaces respectively. Similarly, ξ e ξ v ξ , s i e s i v s i denote the transport of energy. The terms S ξ , S s i are the stress tensors that describe forces acting on a given element due to its surrounding. Therefore the terms S ξ v ξ , S s i v s i yield the mechanical work made by these forces in any macroscopic element of the system. Finally, the terms J ξ q , J s q,i denote the heat fluxes occurring in the system. In case of the bulk phases, where {e k } is an orthonormal basis of R 3 . We do not distinguish the indices below or above since the scalar product of R 3 is given by the metric g kl = δ kl . In (2.12), s ξ kl are functions s ξ kl : D → R, D ⊂ R 3 . On the other hand, S s i are the stress tensors in the tangent bundle of the surfaces Σ i , i.e.
We also assume that all the stress tensors are symmetric (cf. [28,Section 2.2.3]). The fluxes J L q , J G q , J S q ∈ R 3 are the standard vector fields where they are defined and J s q,i (x) ∈ T x Σ i , i = 1, 2, 3, for all x ∈ Σ i . Note that even if we neglect the contribution of the gas density and energy density in the conservation equations, there might be possible fluxes such as J G q . With these definitions, all the terms appearing in equations (2.10) and (2.11) are meaningful in the sense of distributions. The conservation laws (2.7), (2.8), (2.9) are the equations describing the evolution of the full system, after suitable choices of the functions S, J q are made. Before describing these functions, we reformulate the above conservation laws as a set of partial differential equations for the different magnitudes describing the evolution of the system (density, velocity, temperature), both in the bulk and at the interfaces. The conservation laws provide also a set of boundary conditions at the contact line C.

Conservation of mass
Equations (2.7)-(2.9) hold in the sense of distributions and imply corresponding laws in the bulk, on the interfaces and on the contact line. The derivation is rather standard, but in the following we present this derivation for (2.7) in detail for completeness and just state the corresponding results for the other conservation laws.
Let us denote by Ω L , Ω S the liquid and solid domains respectively, i.e. supp χ L ≡ Ω L etc. In the weak formulation, the mass conservation equation (2.7) becomes, for ϕ ∈ C ∞ c (R 3 ×R), (2.13) We calculate II using integration by parts with respect to the space variable. Recalling that S and L are constants, this gives Here we have denoted by n i , i = 1, 2, the unit normal vectors to the surfaces Σ i , i = 1, 2, inward with respect to the liquid domain and by n 3 the unit normal vector to Σ 3 pointing into the gas. We have omitted here the subscript i for the normal vectors for simplicity and will precise wherever it is required. Similarly, the unit vectors ν i are in the plane orthogonal to the contact line C and also tangent to the interfaces Σ i , i = 1, 2, 3, pointing away from the contact line. The surface divergence div Σ is defined in (A.2). Also, we have used the specific form of the solid velocity field (2.1) in the above computation which yields div x v S = 0. Furthermore, we use two transport equations (A.5), (A.6) in the following form, taking into account the choice of the normal vectors, which yield the boundary contributions at the interfaces Σ i (t) and at the moving contact line C(t) respectively, L v s i,n ϕ dt dS, (2.14) similar for the solid domain Ω S (t) and Here, κ Σ denotes the mean curvature of the interfaces and the normal derivative ∂ Σ t is defined in (A.1). In addition, we assume that the normal component of the velocity of the interfaces in contact with the solid is equal to the normal velocity of the solid, i.e.
This means that there is no flux of matter from the solid into the corresponding interfaces.
Returning to the equation (2.13) and integrating by parts with respect to the time variable and using (2.14)-(2.16), we obtain, Therefore, mass conservation in the liquid and at the interfaces Σ i become, (2.19) Note that the mass balance equation gives no condition for the gas in the limit G = 0. Also there is no contribution of the solid density S to the interface equation (2.19) due to the condition (2.16).
Equation (2.19) for i = 2, that is the mass condition at the liquid-solid interface, is a generalization of the standard boundary condition v L n = v s 2,n = v S n for an immiscible fluid. We also obtain the boundary condition for the function s i at the contact line C, derived from the mass balance equation (2.7), or equivalently from (2.17) as, (2.20)

Conservation of momentum
Similarly, we obtain from the momentum balance equation (2.8) the following equations in the bulk and at the interfaces Σ i , which can be obtained from (2.23) by taking s 2 = 0 and v L n = v s 2,n . Further, the momentum balance equation (2.8) implies at the contact line C, Condition (2.25) is equivalent to Young's law in the standard case.

Conservation of energy
Next, the energy equations in the bulk and at the interfaces follow from the energy conservation law (2.9), Recall that e ξ , e s i in the above equations are the total energies. From (2.26), with the help of the balance equation for the kinetic energy in the liquid, which can be obtained by multiplying the momentum equation (2.21) by v L , together with the mass equation (2.18) as (similarly for the solid), one obtains the balance equations for the internal energy u L and u S as, Similarly, subtracting the surface mass equation (2.19) multiplied by 1 2 |v s i | 2 from the surface momentum equation (2.24) multiplied by v s i yields the balance equations for the kinetic energy on the interfaces, which is appearing in the above equation for the surface kinetic energy. We also find the following equation on the liquid-gas interface, Then the following balance equations for the internal energy u s i are obtained from (2.28)-(2.29), with the help of the above kinetic energy equations, and The last term of equation (2.31) appears due to the difference of the velocities of the interface and the liquid which transforms into the internal energy.
The energy equation at the contact line C becomes, from the conservation equation (2.9), The above equation (

Constitutive equations
We now use Onsager's method to derive relations between the generalized fluxes and generalized forces. To this end, we first write the equation for the balance of entropy.

Entropy balance
The general form of the equation satisfied by the specific entropy is, where the total specific entropy is expressed as, given the previous expression for the energy density (i.e. neglecting the contribution of the gas), the total entropy flux as, and the total entropy production as, Here s ξ , s s i are the specific entropies per unit mass, J ξ , J s i are the entropy fluxes per unit mass and σ ξ , σ s i entropy production per unit volume in the bulk and at the interfaces respectively. Note that, even if we neglect the contribution of the gas density and energy density (hence, entropy density) in the conservation equations, there might be possible fluxes such as J G q , J G s . The entropy fluxes are functions of the energy fluxes, heat fluxes etc, discussed previously.
The above balance equation (2.35) yields the classical entropy equations in the bulk, that is, and at the interfaces, and Finally, at the contact line we obtain, The goal here is to obtain the entropy productions σ ξ , σ s i from the energy equations and the Gibbs relation (using the hypothesis of local equilibrium).
Classical thermodynamics proposes the relation dS = dQ T between entropy and heat flux. Combining this with the assumption of local equilibrium, it is natural to assume that In order to compute the production of entropy σ L , σ S , σ G , σ s i , we use the the general Gibbs equation (cf. [10]), where T is the absolute temperature and µ is the chemical potential per unit mass and the thermodynamic pressure p is related to the chemical potential by Taking into account the assumptions that the densities of the liquid and the solid are constants, the above Gibbs relation gives, 44) The meaning of the relations (2.42), (2.43) are precisely that the entropy density of bulk and the interfaces are functions of (internal) energy density and mass density only. Thus we have the constitutive relations for some functions h ξ , h s i and the absolute temperature and the chemical potential in the respective phases are defined as, Here the strategy is to obtain expressions for the entropy productions σ from the entropy equations (2.36)-(2.39) by using the Gibbs relations s = h( , u) in each phase. Writing the convective derivatives (∂ t + v · ∇) (·) with the help of the relation (2.46) and then plugging in the mass and energy equations for , u yields σ in terms of the generalized flux-force pairs, where one further employs Onsager principle.
In the bulk, we therefore write down the entropy production (2.36) with the help of Gibbs relation (2.42) and (2.44) (or, (2.45)-(2.46)) and the flux relation (2.41), using energy equation (2.30) and the incompressibility condition (2.18), as We denote here and below by S ξ := (S ξ + P ξ I), S s i := (S s i + P s i (I − n i ⊗ n i )) the traceless part of the stress tensors where P ξ := − 1 3 tr S ξ , P s i := − 1 2 tr S s i are the mechanical pressures in the bulk and at the interfaces respectively. The factor 1 3 in the mechanical pressure is appearing because of the three-dimensional bulk phases and the factor 1 2 appears due to the two-dimensional interfaces. We also denote by Dv ξ : In the deduction of (2.48), we have used the fact that Dv S = 0 due to the rigid body motion (2.1).
Similarly we obtain from (2.38)-(2.39) together with (2.43)-(2.44), (2.41) combining (2.31) and (2.19) and at Σ 3 , (2.51) In principle, the mechanical pressures P ξ , P s i could be different from the thermodynamical pressures p ξ , p s i . In the case of the liquid, this difference is irrelevant as we assume that it is incompressible. Therefore P L is just Lagrange multiplier that guarantees that div x v L = 0 is satisfied. On the other hand, the density of the interfaces could change and therefore the differences (P s i − p s i ) might be relevant.

Onsager's principle
At this point, we will use Onsager's idea to relate the generalized forces and fluxes which says that, • There is a linear relation between generalized forces and generalized fluxes.
• The generation of entropy must be a non-negative quadratic form.
• The matrix of the coefficients is symmetric in the absence of the magnetic field.
With the help of this principle, one obtains different constitutive equations from the entropy production rate. For example, on the right hand side of (2.47), J L q represents the heat flux and ∇ x 1 T is the force, generated by the gradient of the temperature. Similarly, for the second term, S L denotes the viscous stress flux and Dv L is the force due to the velocity gradient. Hence, Onsager's principle yields in Ω L , Here we further employ the assumption that the effects of the cross-diffusion terms are negligible, i.e. a T,v = a v,T = 0. In the present scenario, we therefore obtain, from (2.47), (2.49), the usual form of the stress tensors in the bulk, and the Fourier's law, Similarly, from the interface equations (2.50)-(2.51), the following relations are deduced, neglecting the cross-diffusion effects: Here, the constants a ξ , b ξ , c ξ , a s i , b s i , c s i , k i , α i , β i are non-negative.

Further constitutive assumptions
In this section, we make some further assumptions that lead to some simplifications in the model.
First of all, we assume that the ratio of the viscosities of the gas and the liquid is of the same order as the ratio of the corresponding densities 2 , that is, Hence, in (2.53), we can assume that b G = 0. Thus, from the momentum equation (2.22), one obtains, under the assumption that c G = 0 (which is the case for an ideal gas with Interface stress tensors The two-dimensional interface can be regarded as hydrodynamically ideal (cf. [28, pp. 195]), i.e. there is no contribution of the viscous forces and thus the stress tensor is given by the pressure only. Therefore, equation (2.56) reduces to, Here we assume that b s i = 0 in (2.56) and c s i = 0 in (2.57), hence P s i = p s i . The surface pressure p s i is defined as the 'surface tension' taken with the negative sign.
Interface mass balance equation In (2.58), the last term on the RHS (hydrodynamic contribution) can be neglected compared to the first term on the RHS (chemical potential) (cf. [28, pp. 196]). Also, in equilibrium, the following must hold, where s ie is the equilibrium surface density. Further, assuming that µ s i can be expanded in a Taylor series around s ie , i.e.
(2.58) can be written as, where τ i is a relaxation time. Note that although the same notation is used to denote the relaxation time and the tangent vectors, we use τ as indices over a vector to indicate tangential component, so there should not be any ambiguity. Using (2.68), we can now simplify the surface mass balance equation (2.19) as can be neglected compared to the capillary effect div Σ S s i (cf. [28, pp. 197]). Also the contribution of the term L v L (v L n − v s i,n ) on the left hand side can be neglected (cf. [28, pp. 198]), leading to the form, which can be splitted into the normal and tangential component, invoking the expression for the liquid stress tensor (2.52), relation (2.65) and the incompressibility assumption of the liquid, as Barotropic condition Lastly, assuming that the interface formation process can be modelled as a barotropic process, a linear equation of state is considered which reflects the changes of the surface pressure due to the compression or extreme rarefaction of the surface phase, The parameter γ i is a constant describing compressibility of the fluid at the interface and s i,0 is the surface density corresponding to zero surface pressure.

Liquid-solid interface
In a similar manner, we can write the boundary conditions on the liquid-solid interface Σ 2 as well, such as (2.69), (2.68), (2.74), together with (2.16). The different condition is a generalized Navier slip condition (instead of (2.73)) which one obtains from (2.59) and (2.70) eliminating S S , together with (2.60). This condition shows that there can be apparent slip on the solid substrate.
Gas-solid interface Also, we assume that the density and the chemical potential of the gas-solid interface do not make any contribution to the liquid and contact line dynamics, therefore, However, the contribution of the surface tension p s 3 remains.
Equations at the contact line Finally, for the equations at the contact line, taking into account the above assumptions, the mass equation (2.20) reduces to, Next we re-write the momentum balance equation (2.25). Let θ i be the angles at which the interfaces Σ i intersect the contact line. Therefore, Then, as the convective momentum fluxes along the interfaces are negligible compared to the surface tensions, we can write (2.25) as, which gives taking projecting on ν 2 direction, In the particular situation when θ 2 = π, the above relation yields the Young equation involving the contact angle θ between the liquid-gas and liquid-solid interface, Further incorporating the interface stress relation (2.64), the above equation reduces to, The square term can be neglected compared to the chemical potential in the above relation, due to the similar argument as used before to deduce (2.68) from (2.58) (cf. [28, pp. 196]). Moreover, employing the assumptions concerning the gas-solid interface (2.76) and the mass balance condition at the triple junction (2.77), relation (2.79) yields, But we also have the following expression for the chemical potential in terms of the density (cf. (2.66), (2.67)), Therefore, we obtain from the jump condition for the chemical potential (2.80), with the assumption that m 2 = 0 (i.e. the chemical potential varies from its equilibrium value), But we can further simplify the above relation with the help of the mass equation (2.77) to, where m = m 1 m 2 and under the assumption that s 1 = 0 and (v s 1,τ − v c ) · ν 1 = 0.

Isothermal limit
We investigate now the situation in the case of an isothermal process more carefully. In particular, the goal is to show that in this regime, the temperature is almost constant in the liquid and consequently, the liquid energy equation and liquid-solid and the liquid-gas surface energy equations do not contribute to the model any more.
First, let us consider that the temperature of the liquid can be expressed as where T 0 is a constant temperature, θ L is a small deviation from the constant temperature, andã L is a non-dimensional quantity given by, .
(2.82) We recall that a L is the heat conductivity of the liquid, appearing in (2.54), L is the characteristic length scale of the underlying system, c L m is the specific heat capacity of the liquid (cf. (2.6)), and [·] denotes the dimension of a physical quantity. The first quantity of the RHS in (2.82) represents the ratio of the effect of the thermal conductivity to the source term due to the viscous forces in the equation (2.30), while the second term is the ratio of the thermal conductivity to the transport of heat due to convection. The other terms in the above expression are chosen in such a way in order to make the contribution of the heat flux from the liquid in the surface energy equation (2.31) much larger, compared to all the other terms (as explained below).
In obtaining ( And importantly, the temperature difference between two bulk phases is small compared to the temperature in the liquid dT T . This last assumption is necessary in order to obtain isothermal limit. Otherwise there might be large heat flux which could make the solid temperature far from constant. With the above definition, if we assume thatã L 1, it means that conduction is much bigger than both the convection and the source due to the viscous force in the internal energy equation of the liquid (2.30). Therefore, all the other terms in (2.30) become negligible compared to the heat flux and in the limitã L → ∞, we obtain a heat equation, On the liquid-solid interface Σ 2 , combining the mass density equation (2.19) with the internal energy equation (2.31) for (i, ξ) = (2, S), we obtain, taking into account the specific forms of the heat fluxes (2.54), (2.55), and also the assumption (2.2), For the liquid-gas interface Σ 1 , we can perform similar analysis as above. At this point, note that although we obtained that the pressure of the gas is constant, it does not necessarily imply that the internal energy (and hence, the temperature) is constant too, due to the ideal gas law. One may only have the following relation, where R G is the ideal gas constant. Plugging in the relations for heat fluxes (2.54), (2.55), the energy equation (2.31) for (i, ξ) = (1, G) transforms to, together with the mass equation Defining the non-dimensional parameter as before, , , and lettingã S → ∞, i.e. assuming that the effect of the conduction is much bigger than the convection, the solid energy equation (2.88) further reduces to a heat equation, Similarly, in the limitã S → ∞, the gas-solid energy equation (2.32) becomes a Neumann condition, since we essentially make the contributions of all the other terms negligible compared to the heat flux from the solid with the above definition (2.89) ofã S , From the above discussion, we obtain that the heat equations (2.83), (2.90) in the bulk together with the Neumann conditions (2.87), (2.85) and (2.91) imply that the temperature in the liquid and the solid is constant T 0 . Thus, the energy equations in the bulk and at the interfaces are completely decoupled from the main contact line problem and do not appear in the final model.

Missing condition:
The energy equation at the contact line (2.81), however, remains as it is, even in the isothermal limit. This third boundary condition has not been taken into account in the original model in [25] due to the assumption of isothermal flows (although mentioned in the book [28, p 210] very briefly, claiming that the quantities such as "surface energy" are artificial concepts in a macroscopic model, and introducing one would result in an error; See also pp. 219). As mentioned in the introduction, the missing condition was first pointed out in [4] and was shown to be necessary for well-posedness of the contact line problem (cf. [7, discussion under "final remarks"]). However, our derived missing condition is somewhat different from the one stated in [4, equation (2.23) or (2.24)], but same as in [7, equation (66)].

Summary of the model
To summarize, we recall below the complete system of the moving contact line, for an incompressible isothermal liquid. Henceforth, we re-name some of the constants, such as b L T → 2µ for the dynamic viscosity of the fluid (which must not be confused with the chemical potential, which in any case is not used explicitly furthermore) and β i T → β i , α i T → α i to be consistent with the standard notations.

Lubrication approximation
The general three-dimensional problem can be considerably simplified in the thin film approximation where the ratio of the characteristic length scales in the vertical and horizontal direction to the film is considered as a small parameter. Below we derive different asymptotic limits depending on different orders of magnitude for the slip length parameter and also for different models.

Thin film model for small slip length
First we consider the simple case of spreading of a viscous, incompressible droplet of liquid over a smooth solid surface where the displaced medium is an inviscid gas, in two dimension.
x y y = h(x, t)

Figure 2: thin film model for droplet
The liquid-gas interface is given by y = h(x, t) which is the free surface and the liquidsolid interface is y = 0. The contact points where these interfaces meet are s i (t) = (Λ i (t), 0). We denote by L the characteristic length scale of the droplet, H the characteristic thickness of the film with H L and ε = H L . We are interested in obtaining the limiting model as ε → 0.
In Shikhmurzaev's model (2.92)-(2.95b), we use the notation v L = (u, v), v s i = (u s i , v s i ), i = 1, 2 and we assume v S = 0, i.e. the fluid droplet is on a plane at rest. We omit the index for the pressure of the liquid and just denote it by p. Let us then consider the following scales, Here The time scale is chosen in such a way in order to make the characteristic time of spreading of the droplet of order one. With the time and length scales fixed, it is then standard to re-scale the velocities accordingly. The pressure is scaled as in the standard thin film approximation, in order not to make the first component of the Navier-Stokes equation trivial (cf. (3.1) 3 ). Also, the scaling for difference of the surface densities is made due to the relation with the surface tensions (2.93g) or, (2.94f).

Reduced model
We show below that equations (2.92)-(2.95c) reduce, to leading order terms in ε, to the following system: in (x, y) ∈ (Λ 1 , Λ 2 ) × (0, h): on y = 0: In the above equations, σ 1e σ 3 and all these non-negative constants are O(1) and dimensionless. The derivation of the above system is sketched below.
To derive (3.1)-(3.4c), let us first compute some related terms in rescaled variables. The free surface, i.e. the fluid-gas interface is now given by y = h(x, t), and the unit normal vector, outward with respect to the gas domain and the corresponding tangent vector are given by 5) and the curvature by We compute the rate of strain tensor Therefore, the normal and the tangential stress on the surface y = h(x, t) are given by, and Dv L · n 1 · τ = ε 3 σ 1e Lµ Similarly, the tangential gradient on y = h(x, t) becomes (cf. (A.12)), (3.9) For the differential terms defined along the surface such as ∂ Σ t , div Σ , we can compute them using the explicit forms in Appendix A.3. Also we will not distinguish below betweenp s 1 and p s 1 , if not necessary. On the liquid-solid surface y = 0, corresponding quantities are simpler. The outward unit normal vector with respect to the solid surface is n 2 = (0, 1), whereas τ = (1, 0), κ Σ = 0.
We also use the linear approximation for the surface pressures (2.93g) and (2.94f) in the following form involving the equilibrium density functions s ie , Further, the following assumptions are made, where the first condition means that the relaxation times τ i are very small compared to the macroscopic time scale and the second condition says the deviation of s i from its equilibrium value s ie , on both interfaces, is also small.
Thus in the limit ε → 0, both the time derivative and the non-linear term vanish compared to the viscous term when Reynolds number satisfies ε 2 Re 1 and one obtains (3.1) 2 . Similarly the second component of the Navier-Stokes equation reduces to (3.1) 3 .
Contact line condition: For the equations at the contact line, the velocity of the contact point can be calculated as v c = ∂ t s i = ε 3 σ 1e µ ∂ t Λ i (t), 0 , i = 1, 2. Also, the contact angle being very small, the normal vectors at the contact point satisfy ν 1 = ν 2 = (1, 0). Therefore, the continuity equation (2.95a) reads as which simplifies to (3.4a), taking into account assumption (3.11). Also the contact angle being small, we approximate cos θ ≈ 1− θ 2 2 . Further, let us assume that p s 3 = σ 3e 1 − ε 2 σ 3 ≥ 0 is a constant. Note that the surface tension on the solid-gas interface works in opposite direction to that of the liquid-solid surface. Therefore, as the surface tensions at equilibrium are balanced, i.e. σ 1e + σ 2e + σ 3e = 0, (3.14) the Young's equation (2.95b) can be written as, using (3.13), Using the equilibrium condition (3.14), the above contact angle condition simplifies to (3.4b).
Finally, the missing condition (2.95c) becomes (3.4c) with standard computation. From now on, we omit the bars for convenience.

Derivation of thin film equation
Next, we would like to derive the lubrication equation for the profile h(x, t) from equations (3.1)-(3.4c). We show below that with g = Indeed, from (3.1) 3 and (3.2d), we can first conclude that p is independent of y. Thus, relation (3.1) 2 gives which yields, due to relations (3.2e) and (3.3c) that Thus, (3.17) becomes, together with (3.2d), (3.18) Now, since due to (3.3d) and (3.3b), we have v| y=0 = 0, one gets from the divergence condition (3.1) 1 , Combining this relation with (3.2a) and (3.2c), we can write the equation for the normal velocity of the free surface in the following conservation form, Therefore, with the help of (3.18), we obtain the thin film equation Finally, the mass balance condition at the triple junction (3.4a) gives, with the help of the equations (3.2f) and (3.3e) for u s i , together with the expression (3.18) for u, (3.20) Note that for the equilibrium surface density in the liquid-solid surface, it holds s 2e > s 1e which means q < 1. Also noting that the velocity of the contact point s i (t) is given by, from the thin film equation (3.15) 1 , we may further simplify the above boundary condition (3.20) to (3.16) 1 .
Observe that in the final thin film model (3.15)-(3.16) we have only two parameters: 1 β which represents the ratio of the slip-length and the film thickness and λ 1 which is of order 1 or small.

Remark 1
In the case when the slip length is small compared to the film height, i.e. 1 β 1, we obtain from (3.15) 1 , which can be seen as the counterpart of the classical thin film equation with no-slip boundary condition

Thin film model for large slip length
If we assume that the slip length is large compared to the film thickness i.e.

Remark 2
The above thin film equation can be compared to the standard one with slip boundary condition,

A solid wedge entering a thin liquid film
In this subsection, we briefly discuss a different thin film model, where a wedge-shaped solid enters a thin liquid film (cf. Figure 3) which leads to movement of the contact line and the formation of a meniscus. We briefly discuss here the thin film model in this situation, both, in the standard approach, i.e. with slip-boundary conditions, and in the approach of Shikhmurzaev. We do not give here the details of the derivation, but rather state the limit models and discuss the main differences. We consider a (horizontal, two-dimensional) thin fluid film of thickness H and a wedgeshaped solid which initially touches the film with an angle δ (cf. Figure 4a) and then enters the liquid such thath(t) = g( t t 0 η with η > 0, whereh(t) denotes the distance between the tip of the wedge and the bottom of the liquid domain. Here t 0 is the time scale that characterises the motion of the solid wedge moving down. The particular case of constant velocity of the wedge corresponds to η = 1, while the case of constant acceleration corresponds to η = 2. As a consequence, the contact point moves and makes an angle (π −θ) between the liquid-solid and the liquid-gas interface. The liquid-gas interface is described by y = h(x, t) and the liquid-solid rigid boundary is given by y =h(t) + δ|x|. The free interface makes an angleθ with the horizontal plane (cf. Figure 4b). Therefore, θ =θ + δ and tanθ = −∂ x h. To obtain the lubrication approximation, we need to assume that the two angles θ and δ are very small and of the same order, which means that the contact angle (π − θ) is close to π in this setting. Note that the thin film approximation would not be valid if δ is large. The contact points are given by s i = (Λ i (t), y).
Finally, we assume that the far field condition h → H as |x| → ∞ holds. Since we are in a symmetric configuration with respect to the y-axis, it is enough to analyse the situation x > 0.

Classical approach
In classical fluid mechanics, conditions at contact lines are considered only at equilibrium, hence the contact angle is prescribed. Starting from the incompressible Navier-Stokes equation in the liquid with a slip boundary condition at the fluid-solid interface, one can derive, after an appropriate rescaling, the following thin film equation where k 1 , and k 2 are rescalings of the angles δ and θ, together with the relation where A can be computed from the inner flow via solving a linear ODE. Relation (3.24), together with the thin film equations (3.21)-(3.23) determines fully the dynamics and the position of the contact point.

Shikhmurzaev's approach
For Shikhmurzaev's approach, we use the same scaling as in Section 3.1, only with the adaptation for the quantitiesh, t 0 , and δ, that is, The difference with the thin film model described in Section 3.1 or 3.2 is that, the equations on the liquid-gas interface (2.93a)-(2.93g) are valid for y = h(x, t), x > Λ(t), while the liquid-solid interface conditions (2.94a)-(2.94f) are valid for y =h(t) + δx, x ∈ (0, Λ(t)), with v S = (0, ∂ th ) (at the upper boundary) and for y = 0 with v S = 0 (at the lower boundary). Using the same scaling as in Section 3.1, we only indicate here the differences from the previous classical case. Essentially, it is enough to calculate the boundary conditions at the liquid-solid interface y =h(t) + δx, x ∈ (0, Λ(t)) together with the contact line conditions (2.95a)-(2.95c).
We obtain as before, equations (3.15) in t > 0, x > Λ. On t > 0 and x = Λ, we obtain the following boundary conditions (cf. (3.16) where (3.16) 4 is replaced by (3.22) 1 ), together with the equation (due to the continuity of velocity u at the contact point) Here A is determined as follows: it solves, Also the surface density ζ 2 satisfies, in the inner region, The symmetry assumption of the configuration further gives u = 0 at x = 0 and thus,

Conclusion and discussion
In this paper, we rederive Shikhmurzaev's model for fluid interfaces. An appealing aspect of the model is that it provides a set of equations describing the variables that characterizes the state of the interfaces. In particular, we study in detail under which conditions it is possible to derive isothermal limits. We also conclude, as it was noticed in some previous works ( [4], [7]) that there is a missing boundary condition at the triple junction in the original formulation of Shikhmurzaev model. We further derive in this paper some thin film approximations that provide simpler PDE models than the original one.
Physical assumptions of Shikhmurzaev's model One of the interesting questions is to understand if Shikhmurzaev's model makes any kind of experimental prediction that differs from the traditional models. It is difficult to get a clear insight into the model with so many parameters, even in selected asymptotic regimes. The real values of different parameters are also hard to obtain and compare. We think that the main physical feature that is contained in Shikhmurzaev's model, but not in the classical ones (including Navier slip boundary conditions, precursors, disjoining pressure etc), is that it contains a characteristic time scale that yields the time in which the interface moves towards its equilibrium value. This time scale is assumed to be very short (10 −6 s). Therefore, the effects associated to that scale would be visible only in very fast phenomena. Shikhmurzaev describes the early stages of situations in which a solid object sinks in a fluid (see the discussion below). In the very early stages of this process, the time scale introduced by Shikhmurzaev might yield a different behaviour than the classical models. However, for longer time scales, there might not be any difference between the behaviour of these different models.
In this direction, a possible experiment that could illustrate the differences between Shikhmurzaev models and classical models is the following (although we did not workout the details further). The idea is to study the dynamics of oscillating bubbles, as they appear in sonoluminiscence experiments. Bubbles in a sound field (with frequencies of at least 100,000 Hz) oscillate and the oscillating bubbles have a changing interface. At those frequencies, the presence of the time scale for the relaxation of the interfaces could interfere with the oscillations. We expect to have different behaviours for the classical models and for Shikhmurzaev models in this setup. Due to the spherical symmetry and the lack of contact lines, it seems possible to analyze the model in detail. The analysis of this model would be particularly relevant in order to assert the validity of Shikhmurzaev's model in situations in which the dynamics of the bubble remains as simple as possible (for instance, region in which the spherical shape of the particle retains, shock waves do not arise etc).
Thin film model From the discussion in Section 3.3.1 and Section 3.3.2 on the thin film model of a solid wedge moving downward into a fluid, with two different approaches, we may conclude the following. The classical model says that the contact angle is created instantaneously (cf. equation (3.22) 2 ), while Shikhmurzaev's approach says that the contact angle is created in a continuous manner (cf. equation (3.25) 2 ). This above model of meniscus formation as a general phenomena (not the thin film approximation) has been mentioned in [29,Fig 6].
following operators are defined as, ∇ x a = ∂a ∂x l l , ∇ x u = ∂u k ∂x l kl , div x u ≡ ∇ x · u = l ∂u l ∂x l , div x A = l ∂a kl ∂x l k .
For every vector field w defined on an interface Σ with the unit normal vector n, we decompose it in the normal and tangential components by means of w = w τ + w n n where w τ = (I − n ⊗ n) · w and w n = w · n.
Note that, according to the above notation, w n is a scalar quantity whereas w τ is a vector.

A.2 Differential geometry
We use ∂ Σ t to denote the normal derivative which is defined as follows: suppose a family of surfaces Σ(t) is parametrized as Σ(t) = {x : x = σ(y, t), y ∈ D ⊂ R 2 } where σ has the property ∂ t σ(y, t) = v s n n, with v s n being the normal velocity of the surface. The existence of such a parametrization for every smoothly evolving family of surfaces Σ(t) is easy to obtain. We then define, for a generic function ϕ, defined over Σ, The surface divergence div Σ is defined, for vector fields u defined on Σ with values in the tangent space T Σ, i.e. u(x) ∈ T Σ(x), for all x ∈ Σ, as div Σ u, ϕ Σ := − Σ u · ∇ τ ϕ, ϕ ∈ C ∞ c (Σ), where ∇ τ ϕ = (I − n ⊗ n) · ∇ x ϕ and ϕ is extended by a constant along the normal in a small neighbourhood of Σ.
In the similar way, the surface divergence div Σ for second order tensors is defined, for For a general smooth interface Σ which is expressed as Σ = {x : x = Γ(y, t), y ∈ D ⊂ R 2 } where Γ is a parametrization, the velocity v s of the interface is given by v s = ∂Γ ∂t (y, t) at each point of the interface. If V n denotes the normal velocity of the interface, then V n = v s · n. (A.3) Note that this relation is independent of the parametrization.
We would like to relate γ with the known parametrization h (which is not necessarily a normal parametrization) and then express ∂ Σ t in terms of h. Let us introduce the following notations, for any generic function w and a fixed point (x, y) on the interface y = h(x, 0), W (x, t) := w(γ(x, t), t) = w(γ 1 (x, t), γ 2 (x, t), t), y = h(x, 0) andw (x, t) := w(x, h(x, t), t).
We can further obtain from (A.11), the following formula for tangential gradient, for a function f on the surface y = h(x, t), wheref (x) := f (x, h(x, t)). One can also obtain the above formula by considering any extension f in R 2 and then computing ∇ f · τ , since the final result does not depend on the choice of extension.

B List of variables
The superscripts ξ, s denote the quantities in the bulk and at the interfaces respectively, for ξ = L, S, G, i = 1, 2, 3.