text:
# References

Wigner’s Little Group, Gauge Transformations and Dimensional Descent

Rabin Banerjee111rabin@boson.bose.res.in and Biswajit Chakraborty222biswajit@boson.bose.res.in

S. N. Bose National Centre for Basic Sciences

JD Block, Sector III, Salt Lake City, Calcutta -700 098, India.

Abstract

We propose a technique called dimensional descent to show that Wigner’s little group for massless particles, which acts as a generator of gauge transformation for usual Maxwell theory, has an identical role even for topologically massive gauge theories. The examples of $`BF`$ theory and Maxwell-Chern-Simons theory are analyzed in details.

Wigner’s little group, which is a subgroup of the Poincare group that leaves the four-momentum invariant, is of fundamental importance for both massive and massless particles. In the latter case this group is $`E(2)`$, which is a semi-direct product of $`T(2)`$ (the group of translations in the plane) and $`SO(2)`$. Each of these subgroups has its own role, but here we are concerned with $`T(2)`$ which acts as a generator of gauge transformation in usual gauge theories, like Maxwell or Kalb-Ramond, where the quanta of excitations are massless. There are however instances where gauge invariance coexists with mass - the topologically massive gauge theories like the $`BF`$ theory or the Maxwell-Chern-Simons(MCS) theory in 2+1 dimensions. Recently we had shown that in the $`BF`$ theory, the gauge invariance is generated by a particular representation of $`T(3)E(3)`$, so that one has to go beyond Wigner’s little group. For the MCS model, however, Wigner’s little group (isomorphic to $`\times 𝒵_2)`$ itself suffices, albeit in a different representation. Is there then a systematic method of obtaining these distinct representations and, if possible, connect them to Wigner’s little group for massless particles, which acts as gauge generator for conventional gauge theories?

The object of this paper is to analyze this and related issues. We show that the specific representations of the translation groups that act as gauge generators of topologically massive gauge theories follow naturally by a dimensional descent of the standard representation of Wigner’s little group in one higher dimension. As a byproduct we show the connection between the helicity quantum numbers of massless particles and the spin quantum number of massive particles in one lower dimension.

Let us first recapitulate certain results for the $`BF`$ theory in 3+1 dimensions, whose dynamics is governed by the following Lagrangian density333Greek indices will always denote 3+1 dimensional space-time. Latin indices from the middle of the alphabet(like $`i,j,k`$) will indicate 4+1 dimensions while those from the beginning(like $`a,b,c`$) will correspond to 2+1 dimensions.

$$=\frac{1}{4}F^{\mu \nu }F_{\mu \nu }+\frac{1}{12}H_{\mu \nu \lambda }H^{\mu \nu \lambda }\frac{m}{6}ϵ_{\mu \nu \lambda \rho }H^{\mu \nu \lambda }A^\rho $$
(1)
where,

$$F_{\mu \nu }=_{[\mu }A_{\nu ]};H_{\mu \nu \lambda }=_{[\mu }B_{\nu \lambda ]};\mu =0,1,2,3.$$
(2)
The polarization vector(tensor) associated with the one(two) form fields $`A_\mu (B_{\nu \lambda })`$ are given by,

$$\epsilon ^\mu =i\left(\begin{array}{c}0\\ a_1\\ a_2\\ a_3\end{array}\right);\epsilon =\{\epsilon ^{\nu \lambda }\}=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& a_3& a_2\\ 0& a_3& 0& a_1\\ 0& a_2& a_1& 0\end{array}\right)$$
(3)
corresponding to the momentum four-vector,

$$p^\mu =(m,0,0,0)^T$$
(4)
for a massive quanta at rest. Notice that the entries in (3) satisfy a duality relation connecting the space components of the two matrices.

It might be mentioned that identical expressions for polarization vector(tensor) are obtained for Proca( massive Kalb-Ramond) theory.444Of course the entries in the different matrices are now independent and are not constrained by the duality relation mentioned earlier which is a consequence of the coupling term in the $`BF`$ model. This is not unexpected if one recalls the equivalence of the $`BF`$ model to such theories.

To begin with the present analysis, note that $`T(3)`$, which acts as the generator of gauge transformation in the $`BF`$ model, is an abelian invariant subgroup of $`E(3)`$. The group $`E(3)`$ is Wigner’s little group for a massless particle in 5 dimensions. This suggests that we shall have to consider 5(= 4+1) dimensional space-time. Now an element of Wigner’s little group in 5 dimensions can be written as

$$W_5(p,q,r;\psi ,\varphi ,\eta )=\left(\begin{array}{ccccc}1+\frac{p^2+q^2+r^2}{2}& p& q& r& \frac{p^2+q^2+r^2}{2}\\ p& & & & p\\ q& & R(\psi ,\varphi ,\eta )& & q\\ r& & & & r\\ \frac{p^2+q^2+r^2}{2}& p& q& r& 1\frac{p^2+q^2+r^2}{2}\end{array}\right)$$
(5)
where $`p,q,r`$ are any real numbers, while $`R(\psi ,\varphi ,\eta )SO(3)`$, with $`(\psi ,\varphi ,\eta )`$ being a triplet of Euler angles. The above result can be derived by following the standard treatment. The corresponding element of the translational group $`T(3)`$ can be trivially obtained by setting $`R(\psi ,\varphi ,\eta )`$ to be the identity matrix and will be denoted by $`W_5(p,q,r)=W_5(p,q,r;\mathrm{𝟎})`$.

Let us now consider free Maxwell theory in 5-dimensions,

$$=\frac{1}{4}F^{ij}F_{ij};i,j=0,1,2,3,4$$
(6)
A free photon will have three independent transverse degrees of freedom. Correspondingly the polarization vector $`\epsilon ^i`$ can be brought to the following form,

$$\epsilon ^i=(0,a_1,a_2,a_3,0)^T$$
(7)
if the photon of energy $`\mathrm{`}\omega `$’ is taken to be propagating in the 4-direction so that the energy momentum 5-vector takes the following form

$$p^i=(\omega ,0,0,0,\omega )^T$$
(8)
Note that this automatically satisfies the ‘Lorentz gauge’ $`\epsilon ^ip_i=0`$. If we now suppress the last rows of the column matrices $`\epsilon ^i`$ (7) and $`p^i`$ (8) we end up with the polarization vector and energy-momentum 4-vector of a Proca quanta of mass $`m=\omega `$ in the rest frame in 3+1 dimensions. This is equivalent to applying the projection operator given by the matrix

$$𝒫=diag(1,1,1,1,0)$$
(9)
and deleting the last rows which now have only null entries. This method of getting results in any dimension by starting from the known results in one dimension higher will be called dimensional descent.

Similarly, to reproduce the polarization tensor of massive Kalb-Ramond(KR) theory in 3+1 dimensions let us consider free massless KR model in 5-dimensions,

$$=\frac{1}{12}H^{ijk}H_{ijk}$$
(10)
Proceeding just as in , the polarization matrix $`\epsilon =\{\epsilon ^{ij}\}`$ can be brought to the following maximally reduced form

$$\epsilon =\{\epsilon ^{ij}\}=\left(\begin{array}{ccccc}0& 0& 0& 0& 0\\ 0& 0& \epsilon ^{12}& \epsilon ^{13}& 0\\ 0& \epsilon ^{12}& 0& \epsilon ^{23}& 0\\ 0& \epsilon ^{13}& \epsilon ^{23}& 0& 0\\ 0& 0& 0& 0& 0\end{array}\right)$$
(11)
Again deleting the last row and column, one gets the polarization tensor$`(\overline{\epsilon })`$ (3) of the massive KR model. This is equivalent to applying the projection operator as $`𝒫\epsilon 𝒫`$. Thus the polarization vector and tensor of the Proca and massive KR models, respectively, have been reproduced.

Now coming to the gauge transformation properties of polarization vector (7) and polarization tensor (11) under the translational subgroup $`T(3)`$, let $`W_5(p,q,r)`$ act on these objects one by one. First, acting on $`\epsilon ^i(\text{7})`$, one gets

$$\epsilon ^i\epsilon ^i=W_5(p,q,r)_{}^{i}{}_{j}{}^{}\epsilon ^j=\epsilon ^i+\delta \epsilon ^i=\epsilon ^i+(pa_1+qa_2+ra_3)\frac{p^i}{\omega }$$
(12)
which is indeed a gauge transformation in (4+1) dimensional Maxwell theory. Applying the projection operator $`𝒫`$ (9) on (12) yields

$$\delta \overline{\epsilon }^\mu =𝒫\delta \epsilon ^i=\frac{1}{\omega }(pa_1+qa_2+ra_3)p^\mu $$
(13)
Here $`\overline{\epsilon }^\mu =(0,a_1,a_2,a_3)^T`$ and, modulus the $`i`$-factor, corresponds to the expression in (3). This is precisely how the polarization vector in $`BF`$ theory transforms under gauge transformation. In fact we can write

$$\delta \overline{\epsilon }^\mu =D_{}^{\mu }{}_{\nu }{}^{}(p,q,r)\overline{\epsilon }^\nu \overline{\epsilon }^\mu $$
(14)
where,

$$D(p,q,r)=\left(\begin{array}{cccc}1& p& q& r\\ 0& 1& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right)$$
(15)
Coming next to the polarization matrix $`\epsilon =\{\epsilon ^{ij}\}`$, its transformation law is given by,

$$\epsilon \epsilon ^{}=W_5(p,q,r)\epsilon W_5^T(p,q,r)=\epsilon +\delta \epsilon $$

where,

$$\delta \epsilon =\{\delta \epsilon ^{ij}\}=\left(\begin{array}{ccccc}0& \alpha _1& \alpha _2& \alpha _3& 0\\ \alpha _1& 0& 0& 0& \alpha _1\\ \alpha _2& 0& 0& 0& \alpha _2\\ \alpha _3& 0& 0& 0& \alpha _3\\ 0& \alpha _1& \alpha _2& \alpha _3& 0\end{array}\right)$$
(16)
where, $`\alpha _1=(q\epsilon ^{12}+r\epsilon ^{13}),\alpha _2=(p\epsilon ^{12}r\epsilon ^{23})),\alpha _3=(p\epsilon ^{13}+q\epsilon ^{23})`$. Again this can be easily recognized as a gauge transformation in (4+1) dimensional KR theory involving massless quanta, as $`\delta \epsilon ^{ij}`$ can be expressed as $`(p^if^j(p)p^jf^i(p))`$ with suitable choice for $`f^i(p)`$, where $`p^i`$ is of the form (8). Now applying the projection operator $`𝒫`$ (9), we get the change in the 3+1 dimensional polarization matrix $`\overline{\epsilon }=\{\epsilon ^{\mu \nu }\}`$, by the formula, $`\delta \overline{\epsilon }=𝒫\delta \epsilon 𝒫^T`$. This simply amounts to a deletion of the last row and column of $`\delta \epsilon `$. The result can be expressed more compactly as

$$\delta \overline{\epsilon }=(D\overline{\epsilon }D^T\overline{\epsilon })$$
(17)
where $`D`$ has already been defined in (15). Again this has the precise form of gauge transformation of the polarization matrix of $`BF`$ model, since it can be cast in the form

$$\delta \overline{\epsilon }_{\mu \nu }=i(p_\mu f_\nu (p)p_\nu f_\mu (p))$$
(18)
for a suitable $`f_\mu (p)`$, where the form of $`p_\mu `$ is now given by (4).

Clearly the generators $`T_1=\frac{D}{p},T_2=\frac{D}{q},T_3=\frac{D}{r}`$ provide a Lie algebra basis for the abelian group $`T(3)`$, as they commute mutually. Note that they also satisfy

$$(T_1)^2=(T_2)^2=(T_3)^2=T_1T_2=T_1T_3=T_2T_3=0$$
(19)
so that

$$D(p,q,r)=e^{pT_1+qT_2+rT_3}=1+pT_1+qT_2+rT_3$$
(20)
The change in the polarisation vector $`\overline{\epsilon }`$ (14) can be expressed as the action of a Lie algebra element,

$$\delta \overline{\epsilon }^\mu =(pT_1+qT_2+rT_3)\overline{\epsilon }^\mu $$
(21)
Besides this $`D`$ also preserves the 4-momentum of a particle at rest. This representation of $`D(p,q,r)`$ was earlier introduced in . However, here we have shown it can be connected to the Wigner’s little group for massless particle in 4+1 dimensions through appropriate projection in the intermediate steps, where the massless particles moving in 4+1 dimensions can be associated with a massive particle at rest in 3+1 dimensions. Similarly the polarization vector and tensor of $`BF`$ theory in 3+1 dimensions can be associated with polarization vector and polarization tensor of free Maxwell and KR theories in 4+1 dimensions.

At this stage, one can ask whether the same features of dimensional descent from higher $`3+1`$ dimensions will go through for the topologically massive Maxwell-Chern-Simons theory in 2+1 dimension. For that let us start with free Maxwell theory in 3+1 dimension. This has two transverse degrees of freedom. Correspondingly the polarization vector $`\epsilon ^\mu `$ takes the following maximally reduced form $`\epsilon ^\mu =(0,a,b,0)^T`$, if the 4-momentum $`p^\mu `$ of the photon of energy ‘$`\omega `$’ and propagating in the 3-direction takes the following form $`p^\mu =(\omega ,0,0,\omega )^T`$. The generator of gauge transformation in this case is $`T(2)`$, which is a subgroup of $`E(2)`$. Now the representation of the 3+1 dimensional Wigner’s little group $`W_4(p,q,\varphi )`$ is just obtained from $`W_5(p,q,r;\psi ,\varphi ,\eta )`$ by replacing $`R(\psi ,\varphi ,\eta )`$ by $`R(\varphi )SO(2)`$, by deleting the columns and rows involving $`r`$ and finally setting $`r=0`$ elsewhere. To obtain a ($`4\times 4`$) representation of $`T(2)`$ one just sets $`\varphi =0`$ in $`W_4(p,q,\varphi )`$.

We are now ready to discuss the dimensional descent from 3+1 to 2+1 dimensions. As a first step, this involves reducing the last rows of $`\epsilon ^\mu `$ and $`p^\mu `$ to zero by using the projection operator $`𝒫=diag(1,1,1,0)`$. The descended objects $`\overline{\epsilon }^a=(0,a_1,a_2)^T`$ ($`a=0,1,2`$) and $`p^a=(\omega ,0,0)^T`$ correspond to the polarization vector and momentum 3-vector of a Proca quanta of mass $`\omega `$ at rest, of the corresponding Proca theory

$$=\frac{1}{4}F^{ab}F_{ab}+\frac{\omega ^2}{2}A^aA_a$$
(22)
in 2+1 dimensions. In order to discuss the gauge transformation properties it is essential to provide a $`3\times 3`$ representation of $`T(2)`$ (denoted by $`\overline{D}(p,q)`$) which amounts to deleting the last row and column of $`D(p,q,r)`$ in (15),

$$\overline{D}(p,q)=\left(\begin{array}{ccc}1& p& q\\ 0& 1& 0\\ 0& 0& 1\end{array}\right)$$
(23)
The corresponding generators are given by,

$$\overline{T}_1=\frac{\overline{D}}{p}=\left(\begin{array}{ccc}0& 1& 0\\ 0& 0& 0\\ 0& 0& 0\end{array}\right);\overline{T}_2=\frac{\overline{D}}{q}=\left(\begin{array}{ccc}0& 0& 1\\ 0& 0& 0\\ 0& 0& 0\end{array}\right)$$
(24)
Just as the Proca theory in 3+1 dimensions maps to the $`BF`$ theory, where gauge transformations were discussed, the Proca theory in 2+1 dimensions is actually a doublet of Maxwell-Chern-Simons theories,

$$=_+_{}$$
(25)
where

$$_\pm =\frac{1}{4}F^{ab}F_{ab}\pm \frac{\theta }{2}ϵ^{abc}A_a_bA_c$$
(26)
with each of $`_+`$ or $`_{}`$ being a topologically massive gauge theory. The mass of the MCS quanta is $`\omega =|\theta |`$, where $`\omega `$ is the parameter entering in (22). We can therefore study the gauge transformation generated in this doublet. The polarization vector for $`_\pm `$, with only one degree of freedom for each of $`_+`$ and $`_{}`$, has been found to be

$$\overline{\epsilon }_\pm ^a=\frac{1}{\sqrt{2}}(0,1,\pm i)^T$$
(27)
while the 3-momentum $`p^a`$ obviously takes the same form as in the Proca model. In analogy with (21), here also one can write,

$$\delta \overline{\epsilon }^a=(p\overline{T}_1+q\overline{T}_2)_{}^{a}{}_{b}{}^{}\overline{\epsilon }^b$$
(28)
where $`\overline{\epsilon }^a=(0,a_1,a_2)^T`$ is the polarization vector for the Proca theory in 2+1 dimensions. Had the Proca theory been a gauge theory, (28) would have represented a gauge transformation, as it can be written as

$$\delta \overline{\epsilon }^a=\frac{pa_1+qa_2}{\omega }p^a$$
(29)
But since Proca theory is not a gauge theory, we can only study the gauge transformation properties of each of the doublet $`_\pm `$ (26) individually. First note that the Proca polarization vector $`\overline{\epsilon }^a`$ is just a linear combination of the two real orthonormal canonical vectors $`\epsilon _1`$ and $`\epsilon _2`$ where,

$$\overline{\epsilon }^a=a_1\epsilon _1+a_2\epsilon _2;\epsilon _1=(0,1,0)^T,\epsilon _2=(0,0,1)^T$$
(30)
Correspondingly the generators $`\overline{T}_1`$ and $`\overline{T}_2`$ (24), form an orthonormal basis as they satisfy $`tr(\overline{T}_a^{}\overline{T}_b)=\delta _{ab}`$. Furthermore,

$$\overline{T}_1\epsilon _1=\overline{T}_2\epsilon _2=(1,0,0)^T=\frac{p^a}{\omega },\overline{T}_1\epsilon _2=\overline{T}_2\epsilon _1=0$$
(31)
On the other hand, the polarization vectors $`\overline{\epsilon }_+^a`$ and $`\overline{\epsilon }_{}^a`$ (27) also provide an orthonormal basis(complex) in the plane as

$$(\overline{\epsilon }_+^a)^{}(\overline{\epsilon }_{}^a)=0;(\overline{\epsilon }_+^a)^{}(\overline{\epsilon }_+^a)=(\overline{\epsilon }_{}^a)^{}(\overline{\epsilon }_{}^a)=1$$
(32)
and can be obtained from the above mentioned canonical ones by appropriate $`SU(2)`$ transformation. This suggests that we consider the following orthonormal basis for the Lie algebra of $`T(2)`$:

$$\overline{T}_\pm =\frac{1}{\sqrt{2}}(\overline{T}_1i\overline{T}_2)=\frac{1}{\sqrt{2}}\left(\begin{array}{ccc}0& 1& i\\ 0& 0& 0\\ 0& 0& 0\end{array}\right)$$
(33)
instead of $`\overline{T}_1`$ and $`\overline{T}_2`$. Note that they also satisfy relations similar to the (1-2) basis,

$$tr(\overline{T}_+^{}\overline{T}_+)=tr(\overline{T}_{}^{}\overline{T}_{})=1;tr(\overline{T}_+^{}\overline{T}_{})=0$$
(34)
One can now easily see that

$$\overline{T}_+\epsilon _+=\overline{T}_{}\epsilon _{}=\frac{p^a}{\omega },\overline{T}_+\epsilon _{}=\overline{T}_{}\epsilon _+=0$$
(35)
analogous to (31). Furthermore,

$$\delta \overline{\epsilon }_\pm ^a=\alpha _\pm \overline{T}_\pm \overline{\epsilon }_\pm ^a=\frac{\alpha _\pm }{\omega }p^a$$
(36)
This indicates that $`\overline{T}_\pm `$ \- the generators of the Lie algebra of $`T(2)`$ in the rotated (complex)basis - generate independent gauge transformations in $`_\pm `$ respectively. One therefore can understand how the appropriate representation of the generator of gauge transformation in the doublet of MCS theory can be obtained from higher 3+1 dimensional Wigner’s group through dimensional descent. A finite gauge transformation is obtained by integrating (36) i.e., exponentiating the corresponding Lie algebra element. This gives two representations of Wigner’s little group for massless particles in $`2+1`$ dimensions, which is isomorphic to $`\times 𝒵_2`$, although here we are just considering the component which is connected to the identity,

$$G_\pm (\alpha _\pm )=e^{\alpha _\pm \overline{T}_\pm }=1+\alpha _\pm \overline{T}_\pm =\left(\begin{array}{ccc}1& \frac{\alpha _\pm }{\sqrt{2}}& i\frac{\alpha _\pm }{\sqrt{2}}\\ 0& 1& 0\\ 0& 0& 1\end{array}\right)$$
(37)
Note that $`G_\pm (\alpha _\pm )`$ generates gauge transformation in the doublet $`_\pm `$, and are related by complex conjugation. This complex conjugation is also a symmetry of the doublet. Such representations were earlier considered in .

Another aspect of this doublet structure is illuminated by observing that, in the 3+1 dimensional case, the photon has two helicity states $`(\pm 1)`$, as follows from topological considerations . By our dimensional descent this corresponds to the two spin states of the Proca quanta in 2+1 dimensions. Each of these states is shared amongst the two components of the doublet(25). This is consistent with the fact that an explicit computation of the spin in the MCS theory$`(_\pm )`$ (25) leads to the same values.

We conclude by observing the unified picture that emerges from this analysis. It might be recalled that Wigner’s little group for massive and massless particles can be regarded as having a common origin. This is because the $`E(2)`$ like little group for massless particle is obtained from the $`SO(3)`$ like little group for massive particles in the infinite momentum/zero mass limit. It is essentially a consequence of the group contraction method connecting these two groups. In this paper we showed the common origin of the gauge generators for usual(i.e., massless) or topologically massive gauge theories. The representations in the latter are obtained by a dimensional descent from one higher dimension of the representation in the former. Following the same ideas it was possible to interpret the spin quantum number of massive particles as the helicity quantum number of massless particles in one higher dimension. As a bonus the structure of the Proca model(in 2+1 dimensions) as a doublet of Maxwell-Chern-Simons theories gets manifested. This is distinctive from the Proca model in 3+1 dimensions which gets mapped to the $`BF`$ model and is devoid of any doublet structure. Finally, the fact that the gauge generators annihilated the physical states was very similar to the Dirac Hamiltonian formulation where the Gauss operators, which act as gauge generators, share the same property. It would be worthwhile to pursue this connection between the Wigner’s little group and the Gauss operators.


format: text

text:
# Interaction of two deformable viscous drops under external temperature gradient

## 1 Introduction

It has long been known that a viscous droplet submerged into another immiscible non-isothermal liquid migrates in the direction of the temperature gradient, $`T`$. This droplet’s migration is due to the temperature-induced surface tension gradient $`_s\left(\gamma \right)`$ at the droplet interface $`S`$ and it is addressed as a thermocapillary migration. Young et al  were first to investigate the mechanism of such migration both theoretically and experimentally. These authors found conditions for the stabilization of a buoyant bubble under gravity force by the thermocapillary force which was directed oppositely. They obtained a theoretical prediction for the migration velocity of a single spherical bubble of radius $`a`$, which is placed in a viscous fluid of viscosity $`\eta _1`$, with an imposed constant temperature gradient $`𝐀`$. In the limit of zero Marangoni number $`\left|\frac{\gamma }{T}\right|\frac{a^2A\rho _1}{\eta _{1}^{}{}_{}{}^{2}}`$ (i.e., with negligible convective transport effects) the droplet velocity $`𝐔`$ in the laboratory coordinate frame is related to the uniform temperature gradient by

$$𝐔=\frac{a𝐀}{\eta _1}\left|\frac{\gamma }{T}\right|\overline{U},\overline{U}=\frac{2}{(2+3\eta )(2+\kappa )}.$$
(1)
Here $`\eta =\eta _2/\eta _1`$ and $`\kappa =\kappa _2/\kappa _1`$ are the ratios of viscosities and thermal conductivities between the internal and the surrounding fluids, respectively, $`\left|\frac{\gamma }{T}\right|`$ is a constant derivative of the surface tension with respect to the local temperature $`T`$ on the interface, and $`\overline{U}`$ is a dimensionless migration velocity. In the case when the droplet and the external fluid have equal viscosity and thermal properties $`\overline{U}=2/150.133`$. A review of the further development of the subject can be found in .

Bratukhin  considered the migration of a deformable viscous droplet at the small Marangoni number. Using perturbation techniques in a spherical coordinate system $`(r,\theta )`$ with origin at the center of the droplet, he found an expression for the deformation of the droplet’s shape $`r=R(\theta )`$:

$`R(\theta )`$ $`=`$ $`1+Ma^2s_2P_2(\theta )+O(Ma^3),`$
$`s_2`$ $`=`$ $`{\displaystyle \frac{3}{4\gamma _0}}\left({\displaystyle \frac{5}{2}}\overline{U}^2{\displaystyle \frac{(1\rho )(1+\eta )}{4+\eta }}+{\displaystyle \frac{Pr\overline{U}}{63\chi }}{\displaystyle \frac{49\chi 7\chi \kappa +18\kappa }{(3+2\kappa )(2+\kappa )}}\right).`$ (2)
Here $`Pr=\nu _1/\chi _1`$ is the Prandtl number, $`\rho =\rho _2/\rho _1`$ and $`\chi =\chi _2/\chi _1`$ are the ratios of density and thermal diffusivity between the internal and the surrounding fluids, respectively, $`\gamma _0=\gamma \rho _1a/\eta _1`$ is the dimensionless surface tension corresponding to the center of the drop, and $`P_2(\theta )`$ is the second order Legendre polynomial.

The first term of $`s_2`$ is related to the inertia term in the Navier Stokes equation, while the second one is related to the convective term in the energy transport equation. It is interesting to note that the expansion of the droplet deformation begins with the term of $`O(Ma^2)`$, and that this term vanishes in the case that combines equal densities of the phases ($`\rho =1`$) and negligible convective heat transport ($`Pr=0`$). It follows from  that the deformations are small also for large values of the reference surface tension $`\gamma _0`$ (small capillary number).

Balasubramaniam and Chai  showed that the Young solution satisfies the full Navier-Stokes equations if the temperature field on the spherical interface is proportional to $`\mathrm{cos}\theta `$, i.e. in the case of negligible convective transport. The perturbation of the spherical form constructed in  under the assumption of small capillary number has the same form as .

In real applications, however, fluid drops are not isolated and their motion is influenced by the presence of neighboring particles. This interaction of drops substantially alters the migration velocities and it may cause drops’ deformations even in the case of negligible convective transport. On the other hand, the second mechanism of the deformations’ suppression (large surface tension) is expected to be valid in the multi body case as well.

Most of the numerous theoretical studies of the interaction of bubbles and droplets in the course of their thermocapillary migration were performed under the assumption of nondeformable drops (zero capillary number). Meyyappan et al  and Feuillbois  studied the Marangoni induced axisymmetric migration of two bubbles and Meyyappan and Subramanian  studied the non-axisymmetric case. Acrivos et al , Satrape  and Wang et al  examined ensembles of more than two bubbles. The remarkable finding of these studies is that a bubble migrates in the presence of another bubbles of the same size with exactly the same velocity as when isolated. Bubbles of different sizes affect the motion of the neighboring bubbles. The influence of the presence of a larger bubble on the motion of a smaller one is more pronounced than its influence on the motion of the larger bubble.

These analyses were extended to the case of liquid drops by Anderson , Keh and Chen  and Wei and Subramanian . In this case, the migration velocity of a drop in the presence of another drops differs from the velocity of an isolated drop even when the drops are of equal size. Similar to the case of interacting bubbles, a larger drop exerts more pronounced influence on the motion of a smaller one. The thermocapillary motion of two non-conducting drops was studied by Loewenberg and Davis .

The literature on the thermocapillary motion of deformable drops is limited in contrast to a somewhat related problem concerning the interaction of drops undergoing a buoyancy driven motion. For the latter case numerous recent studies revealed a rich variety of interaction patterns of deformable drops depending on the Bond number and the initial configuration of the system ( see e.g. Manga and Stone , Zinchenko et al ) and the literature cited). Numerical simulations of an axisymmetric buoyancy-driven interaction of a leading drop and a smaller trailing drop were recently reported by Davis . It was demonstrated that the trailing drop elongates considerably due to the hydrodynamic influence of the leading one. Afterwards, depending on the governing parameters, the drops may either separate and return to spherical shape or the trailing drop may be captured by the leading one, or one of the drops may break up.

For the thermocapillary induced motion, the effect of deformability was studied mostly by a perturbation technique assuming small deformations. Recently, Rother and Davis  applied lubrication approximation to study the effect of a slight deformability of the interfaces on the thermocapillary-driven migration of two drops at close proximity.

Zhou and Davis  presented a recent study in which the drops are free to deform in the course of their thermocapillary interaction. They used a boundary-integral technique to study the thermocapillary interaction of a deformable viscous drop with a larger trailing drop making no a priori assumptions regarding the magnitude of deformations. It was demonstrated that in this case the deformations are small and have a small effect on the drops motion at moderate separation. The influence of deformability becomes significant only when the drops are close together. Zhou and Davis  did not adjust the surface tension to follow the migration of the drops in the external temperature field that would be stationary in the absence of a drop. Hence, the boundary conditions that they used were not compatible with the continuous change of the positions of the drops.

The goal of the present work is to extend the analysis of the thermocapillary interaction of deformable drops to the cases of equal-sized drops and the leading larger drop and to study the influence of the deformations on the relative motion of the droplets. We have applied the boundary integral method for the simulations of the thermocapillary induced motion of two deformable viscous drops in the case of moderate capillary numbers and equal viscosity and thermal properties of the dispersed and continuous phases.

## 2 Problem formulation

Consider two drops embedded in an immiscible viscous fluid with a temperature gradient applied along their line of center. In the absence of gravity the drops migrate in the direction of the applied temperature gradient one beyond another due to thermocapillary forces. The trailing and the leading drop and the continuous fluid are marked by indices $`1`$, $`2`$ and $`3`$, respectively. All the fluids are assumed to be Newtonian and incompressible. The interface tension $`\gamma `$ is assumed to depend linearly on temperature $`\gamma (T)=\gamma (T_0)\frac{\gamma }{T}(TT_0)`$. The temperature field approaches a given linear function $`T=Az+T_0`$ far from the dispersed species. The problem is considered in a laboratory coordinate system, see Fig. 1.

Non-dimensional variables are introduced using the following scaling: the radius $`a_1`$ of the first droplet for the length, $`V^{}=\left|\frac{\gamma }{T}\right|a_1A/\eta _3`$ for the velocity, $`a_1/V^{}`$ for the time, $`\left|\frac{\gamma }{T}\right|A`$ for the pressure and $`Aa_1`$ for the temperature. The problem is governed by the following dimensionless parameters: the Marangoni number, $`Ma`$, the Prandtl number, $`Pr`$, and the capillary number, $`Ca`$, where

$$Ma=\left|\frac{\gamma }{T}\right|\frac{a_{1}^{}{}_{}{}^{2}A\rho }{\eta _{3}^{}{}_{}{}^{2}};Pr=\frac{\nu _3}{\chi _3};Ca=\left|\frac{\gamma }{T}\right|\frac{a_1A}{\gamma _0}.$$
(3)
The other dimensionless parameters are the ratios of the material properties of the fluids, $`\eta =\eta _i/\eta _3,\rho _i/\rho _3,\kappa =\kappa _i/\kappa _3`$ and the geometric parameters $`R=a_2/a_1`$ and $`d=/a_1`$, where $``$ is the distance between the drops along their line of centers.

We shall be interested in the case of small Marangoni numbers and moderate values of the other parameters. Following  one can see that, at the leading order in $`Ma`$, the inertial and convective transport effects are negligible, and the motion inside and outside the droplets is governed by the Stokes equations:

$`p_3`$ $`=`$ $`\mathrm{}𝐯_3,𝐱\mathrm{\Omega }_3,`$
$`p_i`$ $`=`$ $`\eta \mathrm{}𝐯_i,𝐱\mathrm{\Omega }_i,i=1,2,`$
$`𝐯_i`$ $`=`$ $`0,𝐱\mathrm{\Omega }_i,i=1,2,3,`$ (4)
where $`p_i`$ and $`𝐯_i`$ are the pressure and velocity field, respectively. The temperature field inside and outside of droplets satisfies the Laplace equation

$$\mathrm{}T_i=0,𝐱\mathrm{\Omega }_i,i=1,2,3.$$
(5)
The continuous fluid is quiescent far from the drops

$$𝐯_30,\text{as}x\mathrm{}.$$
(6)
At the interfaces the velocity field is continuous,

$$𝐯_0(𝐱)=𝐯_i(𝐱),\text{for}𝐱\mathrm{\Sigma }_i,i=1,2,$$
(7)
and the kinematic boundary conditions require that

$$𝐧_i𝐯_i=𝐕_{ni},\text{for}𝐱\mathrm{\Sigma }_i,i=1,2,$$
(8)
where $`V_{ni}`$ is a velocity of the interface $`\mathrm{\Sigma }_i`$ in the normal direction. The tangential stress balance may be written as

$$𝝉_i(𝚷_3𝚷_i)𝐧_i=\frac{T_i(𝐱)}{\tau },𝐱\mathrm{\Sigma }_i,i=1,2,$$
(9)
and the balance of the normal stresses reads

$$𝐧_i(𝚷_3𝚷_i)𝐧_i=\left(\frac{1}{Ca}T_i(𝐱)\right)H_i,𝐱\mathrm{\Sigma }_i,i=1,2,$$
(10)
Here $`𝚷_3=p_3𝐈+𝐓_3`$ and $`𝚷_i=p_i𝐈+\eta 𝐓_i`$ with $`𝐓_i`$ being the rate of deformation tensor, $`𝝉_i`$ is a unit vector tangential to the $`i^{th}`$ interface, and $`H_i=𝐧_i`$ is the mean curvature on the $`i^{th}`$ interface.

The temperature and the heat flux are continuous

$$T_i=T_3,\kappa \frac{T_i}{n}=\frac{T_3}{n}\text{for}𝐱\mathrm{\Sigma }_i,i=1,2,$$
(11)
and at infinity the temperature tends to a given linear function

$$Tz\text{at}|𝐱|\mathrm{}.$$
(12)
## 3 Boundary-integral formulation. Numerical method

The boundary-value problem outlined in the previous section may also be formulated in an equivalent integral form . Following the classical potential theory  we can obtain a system of boundary-integral equations for the temperature and velocities on the interfaces. When the thermal properties and viscosities of the phases are equal for all $`\mathrm{\Omega }_i`$ phases the boundary integral equations simplify considerably, and the solution to the problem acquires an explicit form. Thus, the temperature and the velocities everywhere are given by

$$T(𝐱)=z,$$
(13)
$`𝐯(𝐱)`$ $`=`$ $`{\displaystyle \frac{1}{8\pi }}{\displaystyle _\mathrm{\Sigma }}𝐉\left[\left({\displaystyle \frac{1}{Ca}}T(𝐲)\right)H(𝐲)𝐧+{\displaystyle \frac{T(𝐲)}{\tau }}𝝉\right]𝑑S(𝐲).`$ (14)
Here $`𝐲\mathrm{\Sigma }=\mathrm{\Sigma }_1+\mathrm{\Sigma }_2`$ and the kernels for the single layer potentials for Stokes flow are

$$𝐉(𝝃)=\frac{𝐈}{\xi }+\frac{𝝃𝝃}{\xi ^3},$$

with $`𝝃=𝐱𝐲`$ and $`\xi =|𝝃|`$. Note that the dimensionless surface tension equals $`1/Caz`$, and the natural restriction $`\gamma 0`$ implies that unequality $`z1/Ca`$ hold for all points on the interfaces during the entire time of the process. Applying this unequality to the initial moment $`t=0`$ provides an upper bound for the capillary number

$$Ca<\frac{1}{2R_1+2R_2+d_0},$$
(15)
where $`d_0`$ is an initial separation between the surfaces of the drops.

In the axisymmetric case, the integrations over the azimuthal angle in  can be performed analytically, and the problem is reduced to the computation of a one-dimensional singular boundary-integral for the interfacial velocity

$$𝐯(𝐱)=_\sigma 𝐁(𝐱,𝐲)𝐟(𝐲)𝑑l(𝐲).$$
(16)
Here $`𝐱=(r_x,z_x)`$ and $`𝐲=(r,z)`$ are vectors on an azimuthal $`(r,z)`$ plane in a cylindrical coordinate system, $`\sigma =\sigma _1+\sigma _2`$ is the union of curves describing the drops’ surfaces in this plane, and $`dl(𝐲)`$ denotes the differential arc length of the corresponding contour.

$$𝐟(𝐲)=\left(\frac{1}{Ca}T(𝐲)\right)H𝐧+\frac{T(𝐲)}{\tau }𝝉=f_n𝐧+f_\tau 𝝉$$

is a traction jump across the interface. The complete expressions for components of the kernel $`B(𝐱,𝐲)`$ are given in the Appendix. Note that the kernel  have a logarithmic singularity at $`𝐲𝐱`$. In the absence of Marangoni effect, when the tangential stresses are continuous ( $`T/\tau =0`$), these singularities are often eliminated making use of a subtraction technique , that is based on a well known equality

$$_{\sigma _i}𝐁(𝐱,𝐲)𝐧(𝐲)𝑑l(𝐲)=0.$$
(17)
When the tangential Marangoni stresses associated with the temperature gradients on the drops surfaces are present, this procedure cannot be simply incorporated into the singularity-subtraction formulation. For the purpose of a study of the droplets migration it is sufficient to compute the normal component of surface velocity $`V_n=𝐧(𝐱)𝐯(𝐱)`$. We employed a subtraction technique  in which the normal component of the velocity on the interfaces can be evaluated as

$$V_n=𝐧(𝐱)𝐯(𝐱)=𝐧(𝐱)_\sigma 𝐁(𝐱,𝐲)𝐟(𝐲)𝑑l(𝐲)$$

$$=_\sigma f_\tau (𝐲)𝐧(𝐱)𝐁(𝐱,𝐲)𝝉(𝐲)𝑑l(𝐲)$$

$$+𝐧(𝐱)_\sigma [f_n(𝐲)f_n(𝐱^{}]𝐁(𝐱,𝐲)𝐧(𝐲)dl(𝐲),$$
(18)
where $`𝐱^{}`$ is an arbitrary chosen point in the bulk or on the interface. The choice $`𝐱^{}=𝐱`$ removes a singularity in the integral. The integrand of the first summand in the right-hand side of  tends to a finite limit at $`𝐲𝐱`$ (see the Appendix).

When the separation between the drops surfaces is small, ’near-singularities’ appear in the integrand if the reference point $`𝐱`$ is located in the gap region. Following Loewenberg and Hinch  we facilitate the calculation of integrals at such ’near singular’ using , where $`𝐱^{}`$ is the image point across the gap region to $`𝐱`$. This approach provides a regularization of the ’near-singularity’ for the regions on the surfaces that are in close proximity.

In our numerical computations, the drop boundaries in the azimuthal plane are discretized using marker points and are approximated by cubic splines. The number of points on the two interfaces is equal. The marker points are distributed uniformly with respect to the arc length (non-uniform distribution was not needed in our simulation because the curvature of the interface remains moderate at all time). The normal and tangent unit vectors at the interfaces and the mean curvatures are computed by taking the derivatives of the position vector with respect to the arc length. When the normal velocities on the interfaces are obtained, the positions of the marker points are advances using Euler or Runge Kutta method. After that the cubic spline approximation of the interfaces is reconstructed and the marker points are redistributed on the new surfaces uniformly with respect to the arc length. The Krasny  filtering procedure was used to smooth the spatial frequency spectrum and to remove noise that was introduced by the computation of high derivatives in the code. All Fourier modes below some ’tolerant number’, which depends on the number of marker points were set to $`0`$. At every step the positions of the drops’ centers of mass were obtained. The velocities of the centers of mass were consequently calculated via numerical differentiation with respect to time. The calculations are stopped when the interfaces of the drops intersect or when the surface tension at any point on the leading drop becomes negative.

## 4 Results and discussion

We consider the axisymmetric drift of two deformable droplets with different size. For simplicity we assume that the droplets and the continuous media have identical viscosity, densities and thermal properties. This ensures that the temperature field is continuous and has continuous derivatives across the interfaces. The droplets are initially spherical and start to move from an initial separation distance. The parameters considered include the capillary number $`Ca`$, the initial droplets’ radii $`R_1`$ and $`R_2`$ and the initial distance $`d_0`$ between surfaces of the droplets. During the motion we observe the evolution of the positions and of the form of the droplets, and the flow patterns inside the droplets and in the external fluid. The relative position was measured by the distance between centers of mass $`L(t)`$ and by the distance between the nearest points on the surfaces of the two drops $`d(t)`$. The deformation of a droplet was characterized by the Taylor deformation parameter

$$D=\frac{2r_{max}\mathrm{\Delta }z}{2r_{max}+\mathrm{\Delta }z},$$
(19)
where $`\mathrm{\Delta }z`$ is the length of the drop along the axis of symmetry, and $`r_{max}`$ is the maximum radial dimension. Positive and negative values of $`D`$ correspond to oblate and prolate shapes, respectively.

### 4.1 Migration velocities. The case of non-deformable drops

We have tested our numerical code by comparing calculated drift velocities with available analytic results. As it was mentioned above, there are analytical results  for the thermocapillary induced motion of non-deformable spherical drops in a uniform external temperature gradient. These results can be applied for the axisymmetric cases of drift of drops with different radii $`R_1`$ and $`R_2`$, and for different separation distances between the centers of mass.

For the test of our code we considered a case of small capillary number, $`Ca=0.001`$, in which the expected deformation of initial spherical shape of the drops during the evolution of the positions is negligibly small. The case of well-separated non-deformable drops was consider by Anderson . He showed that for droplets of equal radii, the velocities are equal and have the form

$$U=\left(1+\frac{1}{L^2}\frac{7}{3}\frac{1}{L^6}\right)|𝐔|.$$
(20)
Keh and Chen  considered also drops at small distance and derived a general formula for the velocities

$`U_1`$ $`=`$ $`M_{11}|𝐔_1|+M_{12}|𝐔_2|`$
$`U_2`$ $`=`$ $`M_{21}|𝐔_1|+M_{22}|𝐔_2|`$ (21)
where $`M_{11},M_{12},M_{21},M_{22}`$ are mobility coefficients that were calculated numerically, and $`|𝐔_1|,|𝐔_2|`$ are droplets’ velocities that each drop would have in the absence of another drop.

The velocities of equal droplets calculated using the formulae ,  and our Boundary-Integral code are depicted in Fig. 2 versus the separation distance. The good agreement in all distance-limits is visible. In the case of large separations our calculations, denoted by solid circles, coincide with Anderson’s asymptotic formula (curve d), and in the limit of infinite distance between the drops we recover the value of the drift velocity $`0.133`$ as in the Young-Bratukhin problem . In the limit of a small separation distance our results are in good agreement with the calculations of Keh and Chen as well (dashed curve). A comparison for the case of drops of different radii is presented in Fig. 3. We observe that the case of the drift velocities of the two non-deformable different droplets is in a good agreement with Keh and Chen predictions.

### 4.2 Deformation of equal-sized drops. Effects of Ca and velocity patterns

Consider now the case of larger capillary numbers, where considerable deformations are expected. The time evolution of the Taylor deformation parameter is shown in Fig. 4 for different values of the capillary number and of the initial separation distance. We observe two qualitatively different scenarios of the deformation development. For sufficiently large initial separations and small $`Ca`$ (e.g. curves d and e), the deformation is small and becomes almost stationary after a short transient period. On the other hand, for the relatively small initial separations and large $`Ca`$ (e.g. curves a and b), the evolution of the deformation after the transient period becomes slower, but remains essentially non-steady. During the initial period, with time scale of $`O(Ca)`$, the deformations occur because the non-uniform surface tension is not compatible with the initial spherical form of the drop. The motion during this period is driven by the surface tension, rather than by its gradient, and the deformation is controlled by the normal component of the interfacial stress. At longer time scales the development of the deformations is induced mostly by the changes in the mutual positions of the drops (characterized by the separation distance). In Fig. 5 the deformation parameter is depicted versus separation distance for equal-size drops and different capillary numbers. The initial separation was kept constant at $`L_0=2.05`$. The leading drop becomes more deformed than the trailing one. This can be explained by the fact that the leading drop is located in a hotter region and, thus, it has a lower surface tension. The dependence of well developed deformations on the capillary number is demonstrated in Fig. 6. Here we plotted the deformation parameter of the more deformed leading drop in Fig. 5 when the separation distance reached $`L=2.04`$. The dependence on small $`Ca`$ appears to be linear.

The deviation of the surfaces from the spherical shape $`\sqrt{r_i^2(\theta )+(z_i(\theta )Z_i)^2}R_i`$ is demonstrated in Fig. 7 for different values of time. One can see considerable deformations in the gap region while the outer regions of the drops remain almost spherical. It is evident from Figs. 4, 5 and 7 that the leading drop deforms to an oblate shape while the trailing one becomes prolate. The mechanism of this deformation is shown in Figure 8 where we plot the velocity field for two equal drops at close proximity for time $`t=10`$. The reference frame moves with the rear point on the axis of the trailing drop as marked in the figure. In this reference frame, the drops are almost at rest except for the gap region. The motion is composed of a flow past an aggregate of two drops and the deformation of the gap region in the opposite direction. The elongation of the trailing drop to a prolate shape and the flattening of the leading drop to the oblate shape are clearly evident.

Similar to the deformation, the velocities of the drops’ centers of mass achieve quasi-steady slowly changing values after a short transient period. For two equal drops at $`Ca=0.2`$ these are plotted in Fig. 2 versus the distance between the surfaces (points on the curves $`a`$ and $`c`$). One can see that, at large separations, the velocities of the two equal-sized drops are almost equal and approach the velocities of the non-deformable drops. In contrast to this, for sufficiently small separation, the velocities of the two deformed drops are different. The leading oblate drop moves slower than the trailing prolate one.

The time evolution of the surface-to-surface distance $`d(t)`$ and center-to center distance $`L(t)`$ for $`Ca=0.2`$ and different initial separations is shown in Fig. 9. We observe two different behaviors depending on the initial separation. At small initial separations, surface to surface distance decreases with time (Figure 9 (e)). For larger separations, surface-to surface distance shows a transient behavior in which it grows initially and decreases after some period (Fig. 9(b)). This initial period of repulsion increases with the initial separation, and for $`d_0=0.5`$ (see Fig 9(c)) the surface-to-surface distance grows during the entire time of the process. To understand this behavior recall that the leading drop deforms into an oblate shape, thus, moving apart from the trailing one. The latter, in turn, deforms into a prolate shape, however, its deformations are lower due to a higher surface tension. As a result, the separation distance increases. The difference in the effective surface tension increases with the separation distance, hence, when the distance is sufficiently small, this effect is overcome by the simultaneous approach of the droplets’ centers of mass (Fig. 9 (d) and (f)) and, after a period of repulsion, the drops surfaces approach each other. The period of initial ’repulsion’ becomes wider for larger initial separations, while when the initial distances are small enough, the difference in the rate of deformation becomes negligible, the initial repulsion period vanishes (see Fig. 9 e) and the drops approach during the entire time of the process. Note that the motion of the droplets towards each other is quite slow, and in the cases that we considered drops that were initially very close together do not progress to a complete collision.

### 4.3 Interaction of Drops of Different Size.

The interaction of equal size drops is governed by the interplay of two effects: relative motion caused by the oblate(prolate) shape of the leading(trailing) drop, and more pronounced deformations of the leading drop in the gap region due to its lower surface tension. These effects take place also for unequal drops, but in this case they are combined with the relative motion caused by the different sizes. Similar to the case of equal drops, after a short transient period the migration velocity of the trailing drop is enchanced and the leading drop motion is retarded compared to the case with no deformations. If the leading drop is smaller than the trailing one, the separation distance decreases with time and the drops come to near contact. For this case, the larger drop remains almost spherical, and the smaller one deforms considerably. The form of the drops in close proximity and the velocity field is shown in Fig. 10 (a) for $`Ca=0.2,R_1=1,R_2=0.5`$. The reference frame moves with the rear point of the trailing drop. The motion is a superposition of a flow past an aggregate of the drops and a motion and deformation of the surface of the leading drop. The larger trailing drop deforms mostly in the near gap region. The deformation pattern in this case agrees qualitatively with the one reported by Zhou and Davis .

The situation is different when the trailing drop is smaller than the leading one. In this case the gap between nondeformable drops increases and the drops separate. The deformability may result in a more complex interaction between such drops if the initial separation is small enough. The velocity field is shown in Fig. 10 (b) for $`Ca=0.2,R_1=1,R_2=1.5`$. The reference frame moves with the point at the center of the gap region. One can observe a very slow motion in the gap and considerable normal velocities at the interfaces in the outer regions, indicating the developing deformations. The evolution of the separation distances and the deformation parameter are shown in Figs 11, 12 and 13 for various drop size ratio and $`Ca`$. Different curves on each plot correspond to different initial separations. It is seen that for the large initial separations, the surface-to-surface distance increases during the entire duration of the process, as in the case when the deformations are absent. When the initial separation is small the distance is first reduced and only increase at the later stages. The transition is controled by the capillary number (see Fig.14) and it seems to occur at separation distance of $`O(10^2)`$.

To understand this unusual behavior, consider the evolution of the deformation factor that is shown in these figures. The smaller drop deforms much stronger than the larger one assuming a prolate form. This is an anticipated result since the deformation of the trailing drop is induced by the pulling presence of the leading one and vice versa. The larger drop induces larger deformations. Since in our case the leading drop is larger, higher deformations of the trailing drop are expected. At the depicted size ratios this effect overcomes the effect of a lower surface tension over the leading drop that has an opposite effect in the case of equal drops, where the leading drop is more deformed.

The deformation of the smaller drop develops very fast during an initial transient period and becomes slower afterwards. The interface of a larger drop is initially pushed by the trailing one and assumes an oblate shape. These deformations are combined with the relative motion of the drops apart from each other due to the different radii. Unless the drops are very close in size, after the initial period the elongation of the trailing drop becomes slower, and the distance between the interfaces begins to grow. The drops separate and the influence of the trailing drop on the leading one causes its deformation to reverse itself and it also becomes prolate (see Figs. 12(c) and 13(c)). When the drops are very close in size (see Fig. 11(c)), the evolution their shape is more symmetric and resembles the one depicted in Fig. 4 for the equal-size case.

During the initial period the drops may come very close together and attractive London van der Waals forces may results in the rupture of the inter-droplet film and cause an eventual coalescence. This process, however, is beyond the scope of our work. The distance between the centers of mass of the drops increases monotonically while the surface exhibits a transient behavior. During the early stages of the dynamics the gap shrinks and the surfaces come to close proximity, while at later stages the drops separate as in the case with no deformations. The time in which the initial perturbation of the shape is relaxed and the gap acquires a minimum appears to depend linearly on the capillary number. In this case the system may result in a higher coalescence rate due to diminishing gap in the transient period.

## 5 Concluding remarks

The interaction between two drops in a temperature gradient shows a variety of behavior patterns. Equal drops always approach each other, however, the surface-to-surface distance has a transient behavior during which the separation first grows and later diminishes. When two unequal drops are moving, with the bigger one leading, an opposite transient behavior is observed if the initial separation distance is small. Here, the surfaces first approach each other and later depart while the drops are moving away from each other. This non-trivial behavior takes place during an initial period of the process. Since our calculations are based on the creeping flow approximation that neglect inertia of the fluid, these results are applicable only if the flow field in the fluid is established much faster than an equilibrium form of the deformed surface. Three different time scales are present in the process under consideration: a viscous relaxation time $`t_\nu =a_1^2\rho /\eta _3`$, capillary surface relaxation time $`t_s=\eta _3a_1/\gamma _0`$, and the thermocapillary time $`t^{}=\eta _3/\left(\frac{\gamma }{T}A\right)`$, that was chosen as the time scale in our considerations. Since $`t_\nu /t^{}=Ma`$ and $`t_s/t^{}=Ca`$, the condition $`Ma<<Ca`$ is sufficient to ensure that the transient period can be described in a creeping flow approximation.

It is interesting to compare the magnitude of deformation in our case with different cases such as buoyancy-driven motion  and that resulting from non-linear dependence of surface tension on the external field . We found that the deformation in the case of a thermocapillary migration causes behavior of the drops at close proximity that is qualitatively different from that observed in the other cases. The calculated magnitude in the thermocapillary was relatively small for both leading and trailing drops and for different aspect ratios. The rich variety of interaction pattern typical for a motion driven buoyancy is not observed. There are several reasons why it is so. The first is a well-known fact that the disturbances of the velocity field, resulting from the thermocapillary-induced motion of a single drop, decay with the distance from the particle as $`1/r^3`$. This is much faster than the $`1/r`$ decay in the case of the gravity-driven motion. Thus, the influence of the presence of a drop on the motion and deformation of the other one is expected to be smaller and to decay more rapidly with the separation distance than the one in the gravity-driven motion case. The other important difference of these two processes is that the Bond number

$$Bo=\frac{\mathrm{\Delta }\rho ga_{1}^{}{}_{}{}^{2}}{\gamma _0},$$

that governs the deformability of the drops undergoing gravity-driven migration, may take any values from $`0`$ to infinity, and the interesting behavior is observed when this value is relatively large, $`Bo>3`$ (see ). In contrast to this, the value of the capillary number that governs the deformability of the drops in the course of their thermocapillary-driven migration is restricted from above if only physically relevant cases of positive interfacial tension are considered. In our calculations the non-vanishing values of $`Ca`$ were typically about $`0.2`$.

Note that using the linear profile of a surface tension – temperature dependence results in restriction on the extent of collective drift of drops. Practically there is some moment of time at which the calculation must stop, because after this time the surface tension becomes zero at the front edge of the leading drop. In the case of temperature dependence this situation implies a phase transition. However, when the surface tension depends on concentration of a solute, zero surface tension cannot be achieved. In this case the dependence of surface tension on concentration ( for some pairs fluid-solute) can be divided to two different regions. In the first region the dependence of tension on concentration is strong and can be described by the linear model considered above. In the second region, when the surface tension is diminished, its change with concentration is very slow and it can be approximated by a constant .

We performed a calculation for the case of equal drops using this simple two-region dependence of the interfacial tension letting the linear dependence down the value $`\gamma /\left(Aa_1\frac{d\gamma }{dT}\right)=10^2`$ and keeping it constant beyond this point. When the drops enter the region of low surface tension the deformation pattern changes and more pronounced distortion of the spherical shapes are obtained. This deformation pattern is accompanied by a clear reduction of the drops’ speed and the entire migration is halted. A demonstration of such behavior is given in Figure 15. More detailed dynamics of such interactions are left for the future study.

## Acknowledgements

This research was supported by The Israel Science Foundation founded by the Academy of Science and Humanities. O. M. L. acknowledges the support of the Israel Ministry for Immigrant Absorption. The authors wish to thank A. M. Leshansky for helpful discussions.

## Appendix: Kernels in Eqn.

If $`r_x0`$, i.e. the point $`𝐱`$ is not on the axis of symmetry

$$B_{zz}=2k\left(\frac{r}{r_x}\right)^{1/2}\left[F(k)+\frac{(zz_x)^2}{(r_x+r)^2+(z_xz)^2}E(k)\right],$$

$$B_{zr}=k\frac{zz_x}{(r_xr)^{1/2}}\left[F(k)\frac{r_x^2r^2+(zz_x)^2}{(r_x+r)^2+(z_xz)^2}E(k)\right],$$

$$B_{rz}=k\frac{zz_x}{r_x}\left(\frac{r}{r_x}\right)^{1/2}\left[F(k)+\frac{r_x^2r^2(zz_x)^2}{(r_x+r)^2+(z_xz)^2}E(k)\right],$$

$$B_{rr}=\frac{k}{r_xr}\left(\frac{r}{r_x}\right)^{1/2}\{[r_x^2+r^2+2(z_xz)^2]F(k)$$

$$\frac{2(zz_x)^4+3(zz_x)^2(r_x^2+r^2)+(r_x^2r^2)^2}{(r_x+r)^2+(z_xz)^2}E(k)\},$$

Here $`F`$ and $`E`$ are complete elliptic functions of the first and second kind respectively,

$$k^2=\frac{4r_xr}{(r_x+r)^2+(z_xz)^2},$$

and $`n_z`$ and $`n_r`$ denote the projection of the normal vector in the azimuthal plane on the axis $`z`$ and $`r`$.

If the point $`𝐱`$ is on the axis of symmetry, $`r_x=0`$, $`k=0`$ and the kernels take the form

$$B_{zz}=2\pi r\frac{r^2+2(zz_x)^2}{[r^2+(zz_x)^2]^{3/2}},B_{rz}=0,$$

$$B_{zr}=2\pi \frac{r^2(zz_x)}{[r^2+(zz_x)^2]^{3/2}},B_{rr}=0,$$

Asymptotic expansions as $`𝐲`$ tends to $`𝐱`$ ($`r_x0`$):

$$B_{zz}=\mathrm{ln}[r^2+(zz_x)^2]+c_{zz}+\mathrm{}B_{rz}=c_{rz}+\mathrm{},$$

$$B_{zr}=c_{zr}+\mathrm{},B_{rr}=\mathrm{ln}[r^2+(zz_x)^2]+c_{rr}+\mathrm{},$$

$$𝐧(𝐱)𝐁(𝐱,𝐲)𝝉(𝐲)=c_{rz}n_r(𝐱)\tau _z(𝐱)+c_{zr}n_z(𝐱)\tau _r(𝐱)+\mathrm{},$$

where $`c_{rr}`$, $`c_{rz}`$, $`c_{zr}`$, and $`c_{zz}`$ are constants. If the curve $`\sigma _i`$ is parametrized by the arc length, the kernels are regular for $`r_x=0`$.

Figure legends

1. Geometric sketch of two viscous drops immersed in an external temperature gradient $`T`$ parallel of their axis of symmetry.
2. Center of mass velocities for equal drops $`R_1=R_2=1`$. The solid line $`(d)`$, dashed line $`(b)`$ and the points on these lines are the results of Anderson, Keh and Chen and our computations for $`Ca=0.001`$, respectively. The points $`(a)`$ and $`(c)`$ correspond to leading and trailing deformable drops at $`Ca=0.2`$.
3. Droplet velocities for $`Ca=0.001`$ and different aspect ratio. In the case $`(a,b)`$ $`R_1=1.0`$ and $`R_2=2.0`$. In the case $`(c,d)`$ $`R_1=1.0`$ and $`R_2=0.5`$. The lines $`b`$ and $`d`$ correspond to the larger drop.
4. The evolution of Taylor deformation parameter defined in Eqn. : $`(a)`$ $`Ca=0.2,d_0=0.1`$; $`(b)`$ $`Ca=0.2,d_0=0.5`$; $`(c)`$ $`Ca=0.1,d_0=0.1`$; $`(d)`$ $`Ca=0.1,d_0=0.5`$; $`(e)`$ $`Ca=0.1,d_0=2`$. Positive and negative values correspond to oblate and prolate shapes, respectively.
5. Deformation parameter for equal-size droplets versus the distance between centers of mass: $`(a)`$ $`Ca=0.2,d_0=0.05`$; $`(b)`$ $`Ca=0.1,d_0=0.05`$; $`(c)`$ $`Ca=0.05,d_0=0.05`$. Positive and negative values correspond to oblate and prolate shapes, respectively.
6. Deformation of the leading droplet versus the capillary number at $`L=2.04`$ for equal-size droplets.
7. The local deformations for trailing $`(a)`$ and leading $`(b)`$ equal-size droplets for the case $`Ca=0.2`$ and $`d_0=0.01`$. 1, 2 and 3 correspond to $`0.5`$, $`3.3`$, and 10 units of time, respectively.
8. Velocity patterns for $`Ca=0.2`$, $`R_1=1.0`$ and $`R_2=1.0`$. The reference frame moves with the speed of the marked point.
9. The evolution of the separation and center-to-senter distances for equal-size drops at $`Ca=0.2`$ as function of the initial separation $`(a)`$. $`(b)`$ Transient dynamics for intermediate initial separation. $`(c,d)`$ Evolution of surface-to-surface separation and center-to-center distance for relatively large initial separation. $`(e,f)`$ Evolution of surface-to-surface separation and center-to-center distance for relatively small initial separation.
10. Velocity patterns for $`(a)`$ $`Ca=0.2`$, $`R_1=1.0`$ and $`R_2=0.5`$; $`(b)`$ $`Ca=0.15`$, $`R_1=1.0`$ and $`R_2=1.5`$. The reference frame moves with the speed of the marked point.
11. Interaction between a leading larger drop and a trailing smaller one at close proximity for $`Ca=0.2,R_1=1.0,R_2=1.1`$: $`(a)`$ evolution of center-to-center distance; $`(b)`$ evolution of surface-to surface separation; $`(c)`$ evolution of Taylor deformation parameter. The numbers 1 to 9 denote the different initial separation.
12. Interaction between a leading larger drop and a trailing smaller one at close proximity for $`Ca=0.15,R_1=1.0,R_2=1.5`$: $`(a)`$ evolution of center-to-center distance; evolution of surface-to surface separation; $`(c)`$ evolution of Taylor deformation parameter. The numbers 1 to 6 denote the different initial separation.
13. Interaction between a leading larger drop and a trailing smaller one at close proximity for $`Ca=0.12,R_1=1.0,R_2=2`$: $`(a)`$ evolution of center-to-center distance; $`(b)`$ evolution of surface-to surface separation; $`(c)`$ evolution of Taylor deformation parameter. The numbers 1 to 7 denote the different initial separation.
14. The relaxation time $`tr`$ with minimum separation distance as function of the capillary number.
15. Deformation of equal-size drops migrating into a region of vanishing surface tension. The upper plot shows the dependence of surface tension on distance. The lower plot demonstrates migration and deformation pattern for various time. $`(a)`$ $`t=0`$; $`(b)`$ $`t=5`$; $`(c)`$ $`t=19`$; $`(d)`$ $`t=25`$.