Terminal thermocapillary migration of a droplet at small Reynolds numbers and large Marangoni numbers

TL;DR

This study investigates the steady-state assumptions in thermocapillary droplet migration at small Reynolds and large Marangoni numbers, revealing that steady solutions do not exist and the process remains unsteady.

Contribution

It provides a detailed analysis of the momentum and energy balances, showing the nonexistence of steady temperature solutions under certain conditions.

Findings

Total momentum is conserved at small Reynolds numbers.

A nonconservative thermal flux prevents steady temperature solutions.

The droplet migration process remains unsteady under studied conditions.

Abstract

In this paper, the overall steady-state momentum and energy balances in the thermocapillary migration of a droplet at small Reynolds numbers and large Marangoni numbers are investigated to confirm the quasi-steady state assumption of the system. The droplet is assumed to have a slight axisymmetric deformation from a sphere shape. It is shown that under the quasi-steady state assumption, the total momentum of the thermocapillary droplet migration system at small Reynolds numbers is conservative. The general solution of the steady momentum equations can be determined with its parameters depending on the temperature fields. However, a nonconservative integral thermal flux across the interface for the steady thermocapillary migration of the droplet at small Reynolds numbers and large Marangoni numbers is identified. The nonconservative integral thermal flux indicates that no solutions of the temperature fields exist for the steady energy equations. The terminal thermocapillary migration of the droplet at small Reynolds numbers and large Marangoni numbers cannot reach a steady state and is thus in an unsteady process.

Equations

$$\begin{array}[]{l}\rho\frac{\partial{\bar{\bf v}}}{\partial t}+\rho\bar{\bf v}\bar{\nabla}\bar{\bf v}=-\bar{\nabla}{\bar{p}}+\mu\bar{\Delta}\bar{\bf v},\\ \gamma\rho\frac{\partial{\bar{\bf v}^{\prime}}}{\partial t}+\gamma\rho\bar{\bf v}^{\prime}\bar{\nabla}\bar{\bf v}^{\prime}=-\bar{\nabla}{\bar{p}}^{\prime}+\alpha\mu\bar{\Delta}\bar{\bf v}^{\prime},\end{array}$$

$$\bar{\mathbf{v}}\to 0,\bar{p}\to p_{\infty },$$

$${\bar{p}}^{\prime }=p_0^{\prime }+P_0^{\prime }(t)$$

$$\begin{array}[]{l}\bar{\bf v}({\bar{\chi}}_{b},{\bar{z}}_{b},t)=\bar{\bf v^{\prime}}({\bar{\chi}}_{b},{\bar{z}}_{b},t),\\ {\bf n}\cdot{\bar{\bf\Pi}}-{\bf n}\cdot{\bar{\bf\Pi^{\prime}}}=\sigma{\bf n}(\bar{\nabla}_{s}\cdot{\bf n})-\bar{\nabla}_{s}\sigma,\end{array}$$

$$\bar{\mathbf{r}}=\mathbf{r}+V_{\infty }t\mathbf{k}$$

$$\begin{array}[]{lll}\bar{\bf v}(\bar{\bf r},t)={\bf v}({\bf r})+V_{\infty}{\bf k},&\bar{p}=p+p_{\infty},&\bar{T}(\bar{\bf r},t)=T({\bf r})+T_{0}+GV_{\infty}t,\\ \bar{\bf v^{\prime}}(\bar{\bf r},t)={\bf v^{\prime}}({\bf r})+V_{\infty}{\bf k},&\bar{p}^{\prime}=p^{\prime}+P^{\prime}_{0}(t),&\bar{T^{\prime}}(\bar{\bf r},t)=T^{\prime}({\bf r})+T_{0}+GV_{\infty}t.\end{array}$$

$$\begin{array}[]{l}\nabla\cdot{\bf vv}=\nabla\cdot{\bf\Pi},\\ \gamma\nabla\cdot{\bf v^{\prime}v^{\prime}}=\nabla\cdot{\bf\Pi^{\prime}},\end{array}$$

$$Re=\frac{\rho v_0{R}_0}{\mu }.$$

$$(v_r,v_{\theta })\to (-V_{\infty }\cos \theta ,V_{\infty }\sin \theta ),p\to 0$$

$$p^{\prime }=p_0^{\prime }$$

$$\begin{array}[]{l}v_{n}(R,\theta)=v^{\prime}_{n}(R,\theta)=0,\\ v_{\tau}(R,\theta)=v^{\prime}_{\tau}(R,\theta),\\ {\bf n}\cdot{\bf\Pi}-{\bf n}\cdot{\bf\Pi^{\prime}}=\frac{1}{We}[\sigma^{*}{\bf n}(\nabla_{s}\cdot{\bf n})-\nabla_{s}\sigma^{*}]\end{array}$$

$$Ca=\frac{v_0\mu }{{\sigma }_0}$$

$$We=ReCa=\frac{\rho v_0^2{R}_0}{{\sigma }_0}.$$

$$\frac{d}{dt}{\int }_V\mathbf{\Xi }dV=\frac{d}{dt}{\int }_W\mathbf{v}dV+\gamma \frac{d}{dt}{\int }_{W^{\prime }}{\mathbf{v}}^{\prime }dV.$$

$$\begin{array}[]{l}\frac{d}{dt}\int_{W}{\bf v}dV=\int_{W}\nabla\cdot{\bf vv}dV=\int_{S_{\infty}}{\bf n}\cdot{\bf\Pi}dS-\int_{S_{1}}{\bf n}\cdot{\bf\Pi}dS,\\ \gamma\frac{d}{dt}\int_{W^{\prime}}{\bf v^{\prime}}dV=\gamma\int_{W^{\prime}}\nabla\cdot{\bf v^{\prime}v^{\prime}}dV=\int_{S_{1}}{\bf n}\cdot{\bf\Pi^{\prime}}dS,\end{array}$$

$$\begin{array}[]{lll}\frac{d}{dt}\int_{V}{\bf{\Xi}}dV&=&\int_{S_{\infty}}{\bf n}\cdot{\bf\Pi}dS-\int_{S_{1}}({\bf n}\cdot{\bf\Pi}-{\bf n}\cdot{\bf\Pi^{\prime}})dS.\end{array}$$

$$\begin{array}[]{ll}\int_{S_{\infty}}{\bf n}\cdot{\bf\Pi}dS&=\int_{S_{\infty}}\Pi_{rr}dS{\bf i_{r}}+\int_{S_{\infty}}\Pi_{r\theta}dS{\bf i_{\theta}}\\ &=\int_{S_{\infty}}(\Pi_{rr}\sin\theta+\Pi_{r\theta}\cos\theta)dS{\bf i_{\chi}}+\int_{S_{\infty}}(\Pi_{rr}\cos\theta-\Pi_{r\theta}\sin\theta)dS{\bf i_{z}},\end{array}$$

$$\begin{array}[]{ll}v_{r}&=v^{(0)}_{r}+Rev^{(1)}_{r}+o(Re)\\ &=-V_{\infty}(1-\frac{1}{r^{3}})\cos\theta-Re\frac{C}{6r^{2}}(1-\frac{1}{r^{2}})(3\cos^{2}\theta-1)+o(Re),\\ v_{\theta}&=v^{(0)}_{\theta}+Rev^{(1)}_{\theta}+o(Re)\\ &=V_{\infty}(1+\frac{1}{2r^{3}})\sin\theta+Re\frac{C}{3r^{4}}\sin\theta\cos\theta+o(Re),\\ p&=\frac{1}{Re}p^{(0)}+p^{(1)}+o(1)\\ &=-\frac{V^{2}_{\infty}}{4r^{6}}+[\frac{V^{2}_{\infty}}{2}(\frac{1}{r^{3}}-\frac{1}{4r^{6}})-\frac{C}{3r^{3}}](3\cos^{2}\theta-1)+o(1),\end{array}$$

$$\begin{array}[]{lll}&&\int_{S_{\infty}}{\bf n}\cdot{\bf\Pi}dS=\int_{0}^{\pi}({\bf n}\cdot{\bf\Pi})|_{r_{\infty}}r^{2}_{\infty}\sin\theta d\theta\\ &=&\{-r^{2}_{\infty}\int_{0}^{\pi}p\sin^{2}\theta d\theta+\frac{r^{2}_{\infty}}{Re}\int_{0}^{\pi}[2\frac{\partial v_{r}}{\partial r}\sin\theta+(\frac{1}{r}\frac{\partial v_{r}}{\partial\theta}+\frac{\partial v_{\theta}}{\partial r}-\frac{v_{\theta}}{r})\cos\theta]\sin\theta d\theta\}{\bf i_{\chi}}\\ &&-\{r^{2}_{\infty}\int_{0}^{\pi}p\cos\theta\sin\theta d\theta-\frac{r^{2}_{\infty}}{Re}\int_{0}^{\pi}[2\frac{\partial v_{r}}{\partial r}\cos\theta-(\frac{1}{r}\frac{\partial v_{r}}{\partial\theta}+\frac{\partial v_{\theta}}{\partial r}-\frac{v_{\theta}}{r})\sin\theta]\sin\theta d\theta\}{\bf i_{z}}\\ &=&\{-r^{2}_{\infty}\int_{0}^{\pi}p^{(1)}\sin^{2}\theta d\theta\\ &&+r^{2}_{\infty}\int_{0}^{\pi}[2\frac{\partial v^{(1)}_{r}}{\partial r}\sin\theta+(\frac{1}{r}\frac{\partial v^{(1)}_{r}}{\partial\theta}+\frac{\partial v^{(1)}_{\theta}}{\partial r}-\frac{v^{(1)}_{\theta}}{r})\cos\theta]\sin\theta d\theta+o(1)\}{\bf i_{\chi}}\\ &&-\{r^{2}_{\infty}\int_{0}^{\pi}p^{(1)}_{1}\cos\theta\sin\theta d\theta\\ &&-r^{2}_{\infty}\int_{0}^{\pi}[2\frac{\partial v^{(1)}_{r}}{\partial r}\cos\theta-(\frac{1}{r}\frac{\partial v^{(1)}_{r}}{\partial\theta}+\frac{\partial v^{(1)}_{\theta}}{\partial r}-\frac{v^{(1)}_{\theta}}{r})\sin\theta]\sin\theta d\theta+o(1)\}{\bf i_{z}}\\ &=&[\frac{\pi}{2r_{\infty}}(\frac{V^{2}_{\infty}}{8}+\frac{C}{r^{2}_{\infty}}+\frac{31V^{2}_{\infty}}{32r^{3}_{\infty}})+o(1)]{\bf i_{\chi}}+o(1){\bf i_{z}}\\ &\sim&\{O(r^{-1}_{\infty})+o(1)\}{\bf i_{\chi}}+o(1){\bf i_{z}}\\ &\to&0.\end{array}$$

$$\begin{array}[]{l}\int_{S_{1}}({\bf n}\cdot{\bf\Pi}-{\bf n}\cdot{\bf\Pi^{\prime}})dS\\ =\frac{1}{We}\int_{S_{1}}[\sigma^{*}{\bf n}(\nabla_{s}\cdot{\bf n})-\nabla_{s}\sigma^{*}]dS\\ =-\frac{1}{We}\oint_{C^{-}}\sigma^{*}dl{\bf\tau}+\frac{1}{We}\oint_{C^{+}}\sigma^{*}dl{\bf\tau}\\ =0,\end{array}$$

$$\frac{d}{dt}{\int }_V\mathbf{\Xi }dV=0.$$

$$\mathbf{n}\cdot \mathbf{\Pi }-\mathbf{n}\cdot {\mathbf{\Pi }}^{\prime }=\frac{1}{We}(2{\sigma }^{\ast }\mathbf{n}-\frac{\partial {\sigma }^{\ast }}{\partial \theta }\tau ),$$

$$\begin{array}[]{lll}\int_{S_{1}}({\bf n}\cdot{\bf\Pi}-{\bf n}\cdot{\bf\Pi^{\prime}})dS&=&\frac{1}{We}\int_{S_{1}}(2\sigma^{*}{\bf n}-\frac{\partial\sigma^{*}}{\partial\theta}{\bf\tau})dS\\ &=&\frac{1}{We}\int_{S_{1}}(2\sigma^{*}\sin\theta-\frac{\partial\sigma^{*}}{\partial\theta}\cos\theta)dS{\bf i_{\chi}}\\ &&+\frac{1}{We}\int_{S_{1}}(2\sigma^{*}\cos\theta+\frac{\partial\sigma^{*}}{\partial\theta}\sin\theta)dS{\bf i_{z}}\\ &=&\frac{1}{We}\int_{0}^{\pi}\sigma^{*}d\theta{\bf i_{\chi}}.\end{array}$$

$$\frac{d}{dt}{\int }_V\mathbf{\Xi }dV=0.$$

$$\begin{array}[]{l}\frac{\partial{\bar{T}}}{\partial t}+\bar{\bf v}\bar{\nabla}\bar{T}=\kappa\bar{\Delta}\bar{T},\\ \frac{\partial{\bar{T^{\prime}}}}{\partial t}+\bar{\bf v^{\prime}}\bar{\nabla}\bar{T^{\prime}}=\lambda\kappa\bar{\Delta}\bar{T^{\prime}},\end{array}$$

$$\bar{T}\to {\bar{T}}_0+G\bar{z},$$

$$\begin{array}[]{l}\bar{T}({\bar{\chi}}_{b},{\bar{z}}_{b},t)=\bar{T}^{\prime}({\bar{\chi}}_{b},{\bar{z}}_{b},t),\\ \frac{\partial{\bar{T}}}{\partial n}({\bar{\chi}}_{b},{\bar{z}}_{b},t)=\beta\frac{\partial{\bar{T}^{\prime}}}{\partial n}({\bar{\chi}}_{b},{\bar{z}}_{b},t).\end{array}$$

$$\begin{array}[]{l}V_{\infty}+\nabla\cdot({\bf v}T)=\frac{1}{Ma}\Delta T,\\ V_{\infty}+\nabla\cdot({\bf v^{\prime}}T^{\prime})=\frac{\lambda}{Ma}\Delta T^{\prime},\end{array}$$

$$Ma=\frac{v_0{R}_0}{\kappa }.$$

$$T\to r\cos \theta$$

Full text

Keywords  Interfacial tension; Thermocapillary migration of a droplet; Quasi-steady state assumption; Microgravity

1 Introduction

In the microgravity environment, a droplet or a bubble placed in a nonisothermal mother liquid migrates in the direction of increasing temperature as a result of surface tension on the interface. This phenomenon is termed the thermocapillary migration of the droplet or the bubble in the fundamental hydrodynamics[1]. Such nonisothermal interfacial flows are very important in a great variety of both macroscopic and microscopic chemical and biological engineering applications[2, 3, 4, 5, 6]. The pioneering study in this area was carried out by Young et al(1959)[7], who gave an analytical prediction of the bubble migration speed in the limit case of zero Reynolds(Re) number and zero Marangoni(Ma) number. Then, Subramanian(1981)[8] introduced a quasi-steady state assumption and extended the YGB results to small Ma numbers. Under the assumption of quasi-steady state, a series of results of theoretical analyses[9, 10] and numerical simulations[11, 12] for migration speeds of bubbles were obtained and are in agreement with those of experimental investigations[13], independent of Ma numbers in the systems.

Although the thermocapillary bubble migration is very well understood, the thermocapillary droplet migration remains a research focus due to its rather complicated behaviors with respect to the transfer of momentum and energy through the interface of two-phase fluids. More attention was paid to the thermocapillary droplet migration in cases of large Re numbers and large Ma numbers. On the one hand, the theoretical analyses and numerical simulations indicate that in the cases of large Re numbers and large Ma numbers, the migration speed of a droplet increases with the increase of Ma number as was reported by Balasubramaniam & Subramanian(2000)[14] and by Ma et al(1999)[15], which is qualitatively different from the corresponding experimental observations[13, 16]. Both the theoretical analyses and numerical simulations are based on the assumptions of the quasi-steady state and the non-deformation of the droplet. On the other hand, Wu & Hu(2007) pointed out an ill-posedness thermal flux boundary condition for the steady thermocapillary droplet migration at large Re numbers and large Ma numbers[17]. Herrmann et al(2008)[18], Wu & Hu(2012)[19] and Samereh et al(2014)[20] investigated the thermocapillary motion of deformable and non-deformable droplets by the numerical methods and indicated that the assumption of the quasi-steady state was not valid in cases of large Re numbers and large Ma numbers. In particular, it was shown[19, 20] that the rise velocities of the droplet against the time/the vertical position in the numerical simulations are in qualitative agreement with those in the experimental investigations[13, 16]. Yin et al(2012)[21] reported a dramatic increase of the computed stable droplet velocities when $Ma>200$ in numerical simulations, which also supports the above investigations. These phenomena were confirmed by Wu & Hu(2012, 2013)[19, 22] and they are related with the fact that a nonconservative integral flux appears across the interface in the steady thermocapillary droplet migration in cases of large Re numbers and large Ma numbers. To avoid the nonconservative integral flux across the interface, it was suggested that a thermal source is added inside the droplet. Under the assumption of the quasi-steady state, an analytical result for the steady thermocapillary migration of the droplet with the source in cases of large Re numbers and large Ma numbers was obtained[23]. However, the physical mechanism of the thermocapillary droplet migration in case of large Ma numbers is still not very well understood.

In this paper, we focus on the thermocapillary migration of a droplet at small Re numbers and large Ma numbers. It is closely related to the heat transfer in heavy oil with the Prandtl number $Pr\sim O(10^{3})$ and the transport processes involving mass transfer[24, 25]. In general, the overall steady-state momentum and energy balances of two phases in a flow domain require that the changes of the momentum and the energy of the domain are equal to the differences between the total momentum and energy entering the domain and those leaving the domain[26]. Firstly, we confirm the overall steady-state momentum balance of the continuous phase and the fluid in the droplet in the thermocapillary droplet migration in the co-moving frame of reference at small Re numbers. Then, by using the asymptotic expansion method, we investigate the continuity of the integral thermal flux across the interface based on the overall steady-state energy balance in the flow domain, and analyze the quasi-steady thermocapillary migration of the droplet in cases of small Re numbers and large Ma numbers. Finally, the qualitative difference between the terminal thermocapillary migration of a droplet at small Ma numbers and that in case of large Ma numbers is discussed, focusing on its intrinsic physical mechanisms.

2 Overall steady-state momentum balance

Consider the thermocapillary migration of a spherical droplet of radius $R_{0}$, density $\gamma\rho$, dynamic viscosity $\alpha\mu$, thermal conductivity $\beta k$, and thermal diffusivity $\lambda\kappa$ in a continuous phase fluid of infinite extent with density $\rho$, dynamic viscosity $\mu$, thermal conductivity $k$, and thermal diffusivity $\kappa$ under a uniform temperature gradient $G$ in cases of small Re numbers and large Ma numbers. The change rate of the interfacial tension between the droplet and the continuous phase fluid as a function of temperature is denoted by $\sigma_{T}(<0)$. Axisymmetric momentum equations for the continuous phase and for the fluid in the droplet in a laboratory cylindrical coordinate system $({\bar{\chi}},{\bar{z}})$ denoted by a bar are written as follows

[TABLE]

where $\bar{\bf v}$ and $\bar{p}$ are the velocity and the pressure, and a prime denotes quantities inside the droplet. The solutions of Eqs. (1) have to satisfy the boundary conditions at infinity

[TABLE]

at the center of the droplet

[TABLE]

and at the interface $({\bar{\chi}}_{b},{\bar{z}}_{b})$ of the two-phase fluids

[TABLE]

where $p_{\infty}$ is the undisturbed pressure of the continuous phase, $\bar{\bf\Pi}$ and $\bar{\bf\Pi}^{\prime}$ are the stress tensors defined in terms of the local fluid velocity and the pressure fields of the two-phase fluids. ${\bf n}$ and ${\bf\tau}$ are the unit vectors normal to the interface and tangent to its meridian, respectively. $\bar{\nabla}_{s}(=\bar{\nabla}-{\bf n}\frac{\bar{\partial}}{\bar{\partial}n})$ is the surface gradient operator. $\sigma[=\sigma_{0}+\sigma_{T}({\bar{T}}-T_{0})]$ is the coefficient of the interface tension decreasing linearly with the temperature ${\bar{T}}$.

Under the quasi-steady state assumption, the velocity fields are fully established at every moment of time. After an initial unstable migration process, the droplet migration may reach a steady state, i.e., migrating with a constant speed $V_{\infty}$. In this case, the momentum equations (1) with the boundary conditions (2)(3)(4) can be established in a coordinate system moving with the droplet velocity $V_{\infty}$ under the coordinate transformation

[TABLE]

and variable transformations

[TABLE]

The droplet is assumed to have a slight axisymmetric deformation from its original sphere shape. The radius of the droplet $R_{0}$, the velocity $v_{0}=-\sigma_{T}GR_{0}/\mu$ and $GR_{0}$ are used as the reference quantities to make the coordinates, the velocity and the temperature dimensionless. The momentum equations (1) combined with the continuous equations can be rewritten in the following dimensionless form in a spherical coordinate system ($r,\theta$) moving with the droplet velocity $V_{\infty}$

[TABLE]

where ${\bf v}=(v_{r},v_{\theta})$, ${\bf\Pi}=-p{\bf I}+\frac{1}{Re}[\nabla{\bf v}+(\nabla{\bf v})^{T}]$ and ${\bf\Pi^{\prime}}=-p^{\prime}{\bf I}+\frac{\alpha}{Re}[\nabla{\bf v^{\prime}}+(\nabla{\bf v^{\prime}})^{T}]$. The Reynolds number is defined as

[TABLE]

The boundary conditions (2)(3)(4) are rewritten in the dimensionless form as

[TABLE]

at infinity,

[TABLE]

at the center of the droplet and

[TABLE]

at the interface $R[=1+f_{1}(\theta)Re+f_{2}(\theta)Re^{2}+O(Re^{3})]$ of the two-phase fluids. The Capillary and Weber numbers are defined, respectively, as

[TABLE]

and

[TABLE]

The non-dimensional surface tension coefficient is now expressed as $\sigma=1-Ca(T+V_{\infty}t)=\sigma^{*}-CaV_{\infty}t$. $P^{\prime}_{0}(t)=p_{\infty}+\frac{1}{Re}V_{\infty}t\nabla_{s}\cdot{\bf n}\approx p_{\infty}+\frac{2}{Re}V_{\infty}t$. Thus, once the droplet migration reaches a steady state, the above problem (1)-(4) in the laboratory coordinate system can be described by the steady momentum equations (7) with the boundary conditions (9)(10)(11) in the coordinate system moving with the droplet velocity. The overall steady-state momentum balance of the two phases in the flow domain in the co-moving frame of reference is also guaranteed.

To confirm the overall steady-state momentum balance of the two phases under the above-mentioned assumption, the entire system $V$ is divided into two subsystems, i.e., the continuous phase domain $W$ $(r\in[R,r_{\infty}],\theta\in[0,\pi])$ and the droplet region $W^{\prime}$ $(r\in[0,R],\theta\in[0,\pi])$. We can write the time rate of the momentum per unit volume ${\bf\Xi}$ in the system $V$ as

[TABLE]

For the righthand side terms in Eq. (14), we apply the generalized transport theorem[26] and the momentum equations (7)

[TABLE]

where $S_{1}$ is the interface of the two-phase fluids. Substituting these expressions into Eq. (14), we have

[TABLE]

For the axisymmetric droplet, the first term on the righthand side of Eq. (16) is rewritten as

[TABLE]

where $\Pi_{rr}=-p+\frac{2}{Re}\frac{\partial v_{r}}{\partial r}$ and $\Pi_{r\theta}=\frac{1}{Re}(\frac{1}{r}\frac{\partial v_{r}}{\partial\theta}+\frac{\partial v_{\theta}}{\partial r}-\frac{v_{\theta}}{r})$. $({\bf i_{r}},{\bf i_{\theta}})$ and $({\bf i_{\chi}},{\bf i_{z}})$ are the unit vectors in the spherical coordinate system $(r,\theta)$ and the cylindrical one $(\chi,z)$, respectively. Since Eq. (9) gives only an approximation of the velocity and the pressure fields at infinity, we have to obtain their asymptotic expansions with respect to the small Re number. In the Appendix, a general solution of the momentum equations (7) for the thermocapillary migration of the droplet is obtained in the form of an expansion in powers of the Re number. The scaled velocity and pressure fields in the continuous phase are described as

[TABLE]

where $V_{\infty}$ and $C$ are constants. Using the above expressions, Eq. (17) is derived as

[TABLE]

Moreover, the second term in the righthand side of Eq. (16) is simplified by means of the boundary conditions (11) and the Stokes theorem as

[TABLE]

where $l$ is the arc length along the closed contour $C^{-}/C^{+}$ (the left/right limit of the equator $C$) on the interface $S_{1}$. Finally, substituting the expressions (19)(20) into Eq.(16), we have

[TABLE]

This implies the overall steady-state momentum balance of two phases in the flow domain in the co-moving frame of reference. Therefore, under the quasi-steady state assumption, the total momentum of the thermocapillary droplet migration system in case of small Re numbers is conservative.

The quasi-steady state assumption not only implies the existence of the steady velocity fields for Eq. (7), but also requires that the net force acting on the droplet by the continuous phase must vanish, i.e., $\int_{S_{1}}{\bf n}\cdot{\bf\Pi}dS=0$. Substituting it into Eq. (20), we have $\int_{S_{1}}{\bf n}\cdot{\bf\Pi^{\prime}}dS=0$. In this case, using the above result and the expression (19), Eq. (15) is reduced to $\frac{d}{dt}\int_{W}{\bf v}dV=\gamma\frac{d}{dt}\int_{W^{\prime}}{\bf v^{\prime}}dV=0$. Thus, under the quasi-steady state assumption, the momentums in two subsystems of the thermocapillary droplet migration system in the case of small Re numbers are also conservative.

Note that the velocity and the pressure fields at infinity at large Re numbers[27], obtained by truncating the expressions (18) at $O(Re)$ and $O(1)$, respectively, can be taken as a special case of small Re numbers. So, in view of the above derivations, the total momentum of the steady thermocapillary droplet migration system at large Re numbers and the momentums in its two subsystems are also conservative[22].

It is meaningful to investigate whether the overall steady-state momentum balance in the thermocapillary droplet migration is valid when the droplet is assumed to be non-deformable and kept in the sphere shape($R=1$). In this case, the stress balance condition at the interface in Eqs. (11) is simplified as

[TABLE]

where $\nabla_{s}\cdot{\bf n}=\nabla\cdot{\bf n}=2$. Using the above boundary condition, the second term in the righthand side of Eq. (16) becomes

[TABLE]

Due to the axisymmetric features of the interface $S_{1}$, the interaction $\frac{1}{We}\int_{0}^{\pi}\sigma^{*}d\theta{\bf i_{\chi}}$ in the cross section (${\bf i_{\chi}},{\bf i_{z}}$) will be balanced with that in the cross section ($-{\bf i_{\chi}},{\bf i_{z}}$). Taking into account of the balance of the interactions in the axisymmetric interface, we have $\int_{S_{1}}({\bf n}\cdot{\bf\Pi}-{\bf n}\cdot{\bf\Pi^{\prime}})dS$ =0. Thus, substituting the above result and the expression (19) into Eq.(16), we have

[TABLE]

This implies that, under the assumptions of the quasi-steady state and the non-deformation of the droplet, the total momentum of the thermocapillary droplet migration system is also conservative.

3 Overall steady-state energy unbalance

The axisymmetric energy equations for the continuous phase and for the fluid in the droplet in a laboratory cylindrical coordinate system $({\bar{\chi}},{\bar{z}})$ denoted by a bar are written as

[TABLE]

where $\bar{T}$ is the temperature. The solutions of Eqs. (25) have to satisfy the boundary condition at infinity

[TABLE]

where ${\bar{T}}_{0}$ is the undisturbed temperature of the continuous phase and the boundary conditions at the interface $({\bar{\chi}}_{b},{\bar{z}}_{b})$ of the two-phase fluids

[TABLE]

Under the quasi-steady state assumption, the temperature fields are fully established at every moment of time. In a similar way, under the coordinate and variable transformations (5)(6), the energy equations (25) combined with the continuous equations can be rewritten in the following dimensionless form in the spherical coordinate system ($r,\theta$) moving with the droplet velocity $V_{\infty}$

[TABLE]

where the Marangoni number is defined as

[TABLE]

The boundary conditions (26)(27) are rewritten in the dimensionless form as

[TABLE]

at infinity and

[TABLE]

at the interface of the two-phase fluids. In the same way, the assumption can be confirmed by the overall steady-state energy balance of the two phases in the flow domain in the co-moving frame of reference. To obtain the overall steady-state energy transport of two phases in the flow domain, Eq. (28) is integrated in the continuous phase domain $W$ and the droplet region $W^{\prime}$. Following a similar derivations as in [22], we have

[TABLE]

and

[TABLE]

Since Eq. (30) only gives an approximation of the temperature field at infinity, we have to obtain an asymptotic expansion of $T$ for the integration of Eq. (32). In the limit case of large Ma numbers, we may apply the perturbation theory to determine the asymptotic behavior of $T$ at $r\gg 1$. In general, the steady migration velocity $V_{\infty}[=\frac{1}{2+3\alpha}\int_{0}^{\pi}T(1,\theta)\sin\theta\cos\theta d\theta>0]$ is a nonlinear function of the Ma number and has three possible functional relationships with the Ma number in the case of the large Ma numbers.

(A) $V_{\infty}$ is an increasing function of the Ma number, $V_{\infty}=V^{[\varphi]}_{\infty}G(Ma)+o[G(Ma)]$, where $V^{[\varphi]}_{\infty}>0$ and $O[G(Ma)]>O[1]$. The function $G(Ma)$ is of the highest order in the nonlinear functional set to show the relationship of $V_{\infty}$ with Ma, such as $Ma^{\xi}$, $Ma^{\eta}\ln Ma$, $(\ln Ma)^{\zeta}$, $\cdots$, $\xi>0$, $\eta>0$ and $\zeta>0$. Taking $\epsilon=1/\sqrt{V_{\infty}Ma}$ as a small parameter, we have $v_{r}=v^{[0]}_{r}+O[\epsilon]$, $v_{\theta}=v^{[0]}_{\theta}+O[\epsilon]$ and $T=T^{[0]}+O[\epsilon]$. The energy equation for the continuous phase in Eq. (28) truncated at the zeroth order $O[1]$ may be written as

[TABLE]

Introducing $T^{[0]}_{*}=T^{[0]}-r\cos\theta$, we have

[TABLE]

In the case of small Re numbers, the velocity fields $(v_{r},v_{\theta})$ have the expansions against Re in Eqs. (18). This leads to $v^{[0]}_{r}=v^{(0)}_{r}+Rev^{(1)}_{r}+o(Re)$ and $v^{[0]}_{\theta}=v^{(0)}_{\theta}+Rev^{(1)}_{\theta}+o(Re)$. Using the coordinate transformation $\xi=r$ and $\eta=\Psi^{(0)}=-\frac{1}{2}\sin^{2}\theta(r^{2}-1/r)$, Eq.(35) can be expressed as

[TABLE]

where $1/v^{[0]}_{r}\approx 1/v^{(0)}_{r}-Rev^{(1)}_{r}/{v^{(0)}_{r}}^{2}+o(Re)$. The analytical result of the outer temperature field for the continuous phase may be obtained as

[TABLE]

where the operator ”+” before the integral is selected to preserve the monotonously increasing trend of $T(r,0)$ with $r(>1)$ for the continuous phase. Using the expressions (18), we have

[TABLE]

Replacing $\theta$ by $\Psi^{(0)}$ in Eq. (38), we have

[TABLE]

where the operator ”$\pm$” in the integral is to be determined by the value of $\theta$ (the operator $"+"/"-"$ corresponds to $\theta\in[0,\pi/2)/[\pi/2,\pi)$). At $r\gg 1$, the result (39) can be expressed as

[TABLE]

where $\Psi^{(0)}\approx-\frac{1}{2}\sin^{2}\theta r^{2}$. Using the above temperature field at the infinity, Eq. (32) may be simplified as

[TABLE]

From Eq. (33) and Eq. (41), we have

[TABLE]

where $R^{3}=1+3f_{2}(\theta)Re^{2}+O(Re^{3})$. Since both $\beta$ and $\lambda$ are positive, we have

[TABLE]

(B) $V_{\infty}$ is a decreasing function of the Ma number or a constant, $V_{\infty}=V^{[0]}_{\infty}+V^{[\varphi]}_{\infty}G(Ma)+o[G(Ma)]$, where $V^{[0]}_{\infty}>0$ and $O[G(Ma)]<O[1]$. $V^{[0]}_{\infty}$ is related with the zeroth order inner temperature near the interface within the thermal boundary layer, which is to be matched with the zeroth order outer temperature outside the thermal boundary layer[10, 23]. The function $G(Ma)$ is of the lowest order in the nonlinear functional set to show the relationship of $V_{\infty}$ with Ma, such as $Ma^{-\xi}$, $Ma^{-\eta}lnMa$, $(\ln Ma)^{-\zeta}$, $\cdots$, $\xi>0$, $\eta>0$ and $\zeta>0$. Taking $\epsilon=1/\sqrt{Ma}$ as a small parameter, we have $v_{r}=v^{[0]}_{r}+O[\epsilon]$, $v_{\theta}=v^{[0]}_{\theta}+O[\epsilon]$, $T=T^{[0]}+O[\epsilon]$ and $V_{\infty}=V^{[0]}_{\infty}+V^{[\varphi]}_{\infty}G(1/\epsilon^{2})+o[G(1/\epsilon^{2})]$. The energy equation for the continuous phase in Eq. (28) truncated at the zeroth order $O[1]$ may be written as

[TABLE]

In a similar manner, the outer temperature field for the continuous fluid at $r\gg 1$ can be expressed as

[TABLE]

Using the above temperature field at the infinity, Eq. (32) is simplified as

[TABLE]

From Eq. (33) and Eq. (46), we have

[TABLE]

Since both $\beta$ and $\lambda$ are positive, we have

[TABLE]

To sum up, Eq. (43) and Eq. (48) may be rewritten as

[TABLE]

From Eq. (31), the integral thermal flux boundary condition across the interface can be expressed as

[TABLE]

So, if the overall steady-state energy of two phases in the flow domain is balanced, Eq. (49) should be reduced to Eq. (50), which seems impossible. It is termed as a nonconservative integral thermal flux across the interface for the steady thermocapillary droplet migration at small Re numbers and large Ma numbers. This implies the overall steady-state energy unbalance of two phases in the flow domain in the co-moving frame of reference and indicates that the steady temperature fields are impossible. The terminal thermocapillary droplet migration at small Re numbers and large Ma numbers cannot reach a steady state. Its terminal state is thus unsteady.

4 Conclusions and discussions

In summary, to verify the quasi-steady state assumption in the thermocapillary migration of a droplet in cases of small Re numbers and large Ma numbers, we have investigated the overall steady-state momentum and energy balances of the system. The droplet is assumed to have a slight axisymmetric deformation from a sphere shape. Under the quasi-steady state assumption, the total momentum of the thermocapillary droplet migration system at small Re numbers has been shown as conservative. Although the general solution of the steady momentum equations can be determined independently, its parameters still depend on the temperature fields. However, a nonconservative integral thermal flux across the interface for the steady thermocapillary migration of the droplet at small Re numbers and large Ma numbers has been identified. The nonconservative integral thermal flux indicates that there is no solution of the temperature fields for the steady energy equations. The terminal thermocapillary migration of the droplet at small Re numbers and large Ma numbers cannot reach a steady state and is thus in an unsteady process.

From the above analysis and the previous investigations of the thermocapillary droplet migration at large Re numbers and large Ma numbers[22], we can conclude that the terminal thermocapillary migration of a droplet at large Ma numbers is in an unsteady state for the system at any Re numbers. Moreover, it was demonstrated theoretically that the terminal thermocapillary migration of a droplet at small Ma numbers is in a steady state[22]. The qualitative difference between the terminal thermocapillary migration of a droplet at small Ma numbers and that at large Ma numbers is due to the evolution of the mode and the intensity of the heat transfer in the system. In general, when the droplet migrates upward, the thermal convection and conduction are two possible ways of the heat transfer in the continuous phase fluid. Either of them can make a redistribution of the thermal energy around the droplet. However, the droplet can only obtain the thermal energy though the thermal conduction across the interface. The thermal energy is not only transferred into the droplet from the top surface but also out the droplet from the bottom surface. The thermal convection in the droplet can only make a redistribution of the thermal energy. In case of small Ma numbers, the heat transfer in the system duo to the thermal conduction across/around the droplet is stronger than that due to the thermal convection around the droplet, so that both the continuous phase fluid and the droplet can obtain simultaneously the thermal energy. The external and internal temperature of the droplet may have a fixed difference as shown in Eq. (6). This leads to a steady migration process. The steady thermocapillary migration of a droplet in case of small Ma numbers was verified experimentally[7, 28]. At large Ma numbers, the heat transfer in the system due to the thermal convection around the droplet is stronger than that due to the thermal conduction across/around the droplet. During the upward migration process, the increment of the internal temperature of the droplet is far smaller than that of the external temperature of the droplet, so that the difference between the external and internal temperatures of the droplet increases more and more, and cannot reach a fixed value as shown in Eq. (6). This reveals an unsteady migration process. The unsteady thermocapillary migration of a droplet at large Ma numbers was verified experimentally[13, 16]. These intrinsic physical mechanisms of the thermocapillary droplet migration at small numbers and large Ma numbers are conducible to the understanding of the interaction of the droplets in both numerical simulations[29, 30] and experimental investigations [31].

Acknowledgments This research is supported by the National Natural Science Foundation of China through the Grants No. 11172310 and No. 11472284. The author thanks the IMECH research computing facility for assisting in the computation.

Appendix: A general solution of steady momentum equations for thermocapillary migration of a droplet at small Re numbers

The momentum equations (7) of the continuous phase and the fluid in the droplet for the steady thermocapillary migration of a droplet at small Re numbers are rewritten in terms of the stream functions $\Psi(r,\theta)$ and $\Psi^{\prime}(r,\theta)$ as

[TABLE]

where $D^{2}=\frac{\partial}{\partial r^{2}}+\frac{\sin\theta}{r^{2}}\frac{\partial}{\partial\theta}(\frac{1}{\sin\theta}\frac{\partial}{\partial\theta})$. The boundary conditions at infinity and at the center of the droplet can be represented in terms of the stream functions $\Psi(r,\theta)$ and $\Psi^{\prime}(r,\theta)$ as

[TABLE]

In the same way, the boundary conditions (11) at the interface are rewritten as

[TABLE]

where $(v_{r},v_{\theta})=(\frac{1}{r^{2}\sin\theta}\frac{\partial\Psi}{\partial\theta},-\frac{1}{r\sin\theta}\frac{\partial\Psi}{\partial r})$ , $(v^{\prime}_{r},v^{\prime}_{\theta})=(\frac{1}{r^{2}\sin\theta}\frac{\partial\Psi^{\prime}}{\partial\theta},-\frac{1}{r\sin\theta}\frac{\partial\Psi^{\prime}}{\partial r})$ and $\phi(={\rm arctg}\frac{dR/d\theta}{R})$ is an angle between ${\bf n}$ and ${\bf i_{r}}$. Under the quasi-steady state assumption, the steady thermocapillary droplet migration requires that the total net force acting on the droplet is zero. In particular, the zero net force in the vertical direction is expressed as

[TABLE]

In the momentum equations (51) with the boundary conditions (52)-(54), the controlled parameter is the Re number. So the solutions of Eqs. (51) are sought in the forms of expansions in powers of the small parameter Re

[TABLE]

Similarly, the pressure fields in Eqs. (7) have the perturbation expansions in the form

[TABLE]

with the boundary conditions $p\rightarrow 0$ at infinity and $p^{\prime}=p^{\prime}_{0}=\frac{1}{Re}p^{\prime(0)}_{0}+p^{\prime(1)}_{0}+o(1)$ at the center of the droplet. Substituting (55) and (56) into Eqs. (51)-(54), the governing equations and the boundary conditions for the zeroth order solutions $\Psi^{(0)}(r,\theta)$ and $\Psi^{\prime(0)}(r,\theta)$ are, respectively, obtained as

[TABLE]

and

[TABLE]

This is the classical model of the thermocapillary migration of a spherical droplet proposed by Young et al(1959)[7]. The solutions of Eqs. (57) can be determined[32, 33] as

[TABLE]

where $p^{\prime(0)}_{0}$ is a constant, $p^{(0)}(r,\theta)$ and $p^{\prime(0)}(r,\theta)$ are obtained by integrating the zeroth order equations of Eqs. (7) with the boundary condition $p^{(0)}(r_{\infty},\theta)\rightarrow 0$. Meanwhile, the governing equations with the boundary conditions for the first order solutions $\Psi^{(1)}(r,\theta)$ and $\Psi^{\prime(1)}(r,\theta)$ in Eqs. (51) are written as

[TABLE]

and

[TABLE]

where $D^{2}\Psi^{(0)}=0$ and $D^{2}\Psi^{\prime(0)}=-\frac{15}{2}V_{\infty}r^{2}\cos^{2}\theta$. Using the expressions (59), the third condition at the interface in Eqs. (61) is expressed as

[TABLE]

which leads to $\Psi^{(1)}(r,\theta)\sim f_{1}(\theta)\sin^{2}\theta$. The function $f_{1}(\theta)$ is taken as $\omega\cos\theta$ ($\omega$ is a constant) in order to have $D^{2}\Psi^{(1)}\sim f_{1}(\theta)\sin^{2}\theta$. Substituting a trial solution $\Psi^{(1)}(r,\theta)=F(r)\sin^{2}\theta\cos\theta$ for the continuous phase into the first equation in Eqs. (60), we have

[TABLE]

Then, the first condition at infinity in Eq. (61) requires that

[TABLE]

And the third condition at the interface in Eqs. (61), i.e., Eq. (62), demands that

[TABLE]

Similarly, a trial solution of the form

[TABLE]

for the fluid in the droplet is suggested in terms of the fourth condition at the interface in Eqs. (61). Then, the second condition at the center of the droplet in Eqs. (61) requires that

[TABLE]

And the fourth condition at the interface in Eqs. (61) demands that

[TABLE]

Using the expressions (63)(64) and the expressions (66)(67), the fifth condition at the interface in Eqs. (61) is reduced to

[TABLE]

Combining Eq. (68) and Eq. (69), we have

[TABLE]

Furthermore, the first order pressure fields for the continuous phase and for the fluid in the droplet are determined by integrating the first order equations of Eqs. (7) with the boundary condition $p^{(1)}(r_{\infty},\theta)\rightarrow 0$

[TABLE]

where $p^{\prime}_{01}$ is a constant. Using the above stream functions $\Psi^{(0)}$, $\Psi^{(1)}$ and pressure fields $p^{(0)}$, $p^{(1)}$ for the continuous phase, the sixth condition in Eqs. (61) is found to be self-perpetuating.

In addition to the above boundary conditions (53)(54) at the interface, the deformed droplet in the steady thermopicallary migration has to satisfy other conditions based on the following facts. On the one hand, the volume of the deformed droplet remains unchanged, i.e.,

[TABLE]

The condition is also self-perpetuating for the above choice of $f_{1}(\theta)=\omega\cos\theta$. On the other hand, the center of mass of the droplet is always fixed at the origin of coordinates, such as $\int_{V}zdV=0$ in the vertical direction. It is further shown that

[TABLE]

which demands that $\omega=0$, i.e., $f_{1}(\theta)=0$. Finally, the first order stream functions for the continuous phase and for the fluid in the droplet are shown to be

[TABLE]

The above derivations reveal that the functions in the solution of the steady momentum equations at small Re numbers are independent of the temperature fields, but the parameters $V_{\infty}$, $C$, $p^{\prime(0)}_{0}$ and $p^{\prime(1)}_{0}$ in the solution still depend on the temperature fields.