Table Of Content

Dynamic van der Waals Theory of two-phase fluids in heat flow Akira Onuki Department of Physics, Kyoto University, Kyoto 606-8502, Japan (February 2, 2008) We present a dynamic van der Waals theory. It is useful to study phase separation when the temperaturevariesinspace. Weshowthatifheatflowisappliedtoliquidsuspendingagasdroplet 5 at zero gravity, a convective flow occurs such that the temperature gradient within the droplet 0 0 nearlyvanishes. Astheheatfluxisincreased,thedropletbecomesattachedtotheheatedwallthat 2 iswettedbyliquidinequilibrium. Inonecasecorrespondingtopartialwettingbygas,anapparent contact angle can be defined. In the other case with larger heat flux, the droplet completely wets n the heated wall expelling liquid. a J PACS numbers: 05.70.Ln, 64.70.Fx, 44.35.+c, 47.20.Dr 7 ] t In usual theories of phase transitions, the fluctuations ThisrelationfollowsfromthevanderWaalstheory. Our f o of the temperature T are assumed to be small and are schemestartswithEqs.(1),(2),and(4),wherethetemer- s neglected. However, there can be situations in which ature is not yet introduced. We then define the local . t phase transitions occur in inhomogeneous T. For exam- temperature T =T(n,e) by a ple, wetting properties near the gas-liquid critical point m are very sensitive to applied heat flux [1,2] and boiling 1 δ ∂s - = =n , (5) d processes remain largely unexplored [3,4]. To treat such T (cid:18)δeS(cid:19)n (cid:18)∂e(cid:19)n n problemsweproposetostartwithcoarse-grainedentropy o and energy, from which the temperature is defined as a where n is fixed in the derivatives. This definition of T c is analogous to that in a micro-canonical ensemble. We functionaloftheorderparameterandtheenergydensity. [ Let us consider one-component fluids, where the or- also define the local chemical potential µˆ (per particle) 2 der parameter is the number density n. In the van der including the gradient contributions by v Waals theory the energy density e and the entropy s 7 δ M 9 per particle are given by e = dnkBT/2 − ǫv0n2 and µˆ =−T(cid:18)δSn(cid:19) =µ−T∇· T ∇n, (6) 0 s = (d/2)lnT +ln(1/v0n 1)+const., respectively, in eˆ − 1 terms of T andn [5]. Here v0 andǫ representthe molec- where eîn Eq.(3) is fixed, µ= T(∂(ns)/∂n) , and 0 ular volume and the magnitude of the attractive poten- − e 5 tial, respectively, and d is the space dimensionality. To 0 M =CT +K. (7) describesituationsinthepresenceofinterfaces,wegener- / t a alize this descriptionby including gradientcontributions Maximization of under a fixed total particle number m to the total entropy and the total internal energy as drnandafixedStotalenergy leadstothe equilibrium E - Rconditions T =const. and µˆ =const. As first derived by d 1 con S =Z dr(cid:20)ns− 2C|∇n|2(cid:21), (1) vfiwalienthndte=hreWn(saxua)rlfsiasc[5de,e6tt]ee,nrtmshioeinneeqdγubi=lyibMµri(unm,−∞T∞in)−tdexrM(fdadnc2e/ndd/xedn)x2s,2itwy=hpµerrocxe- v: µcx isthechemicalpotentialonthRecoexistencecurve[7]. 1 i = dr e+ K n2 . (2) Next we set up the dynamic equations [8]. Hereafter X E Z (cid:20) 2 |∇ | (cid:21) gravity will be neglected. The mass density ρ = mn (m r being the particlemass)obeys∂ρ/∂t= (ρv), where a The gradient terms represent a decrease of the entropy v is the velocity field. The momentum d−e∇ns·ity ρv obeys andanincreaseofthe energydue to inhomogeneityofn. Here the internal energy density is written as ∂ ρv = (ρvv) (↔Π ↔σ). (8) ∂t −∇· −∇· − 2 eˆ=e+K n /2, (3) |∇ | ↔ The Π = Π is the reversible stress tensor invariant ij including the gradient part. We assume that the coeffi- with respe{ct to}time reversal v v. The ↔σ = σ is ij →− { } cientsC andK aresimplyconstantsandthatsinEq.(1) the dissipative stress tensor of the form σ = η( v + ij i j depends on n and e as v )+ (ζ 2η/d)δ v with = ∂/∂x in∇terms j i ij i i ∇ − ∇ · ∇ of the viscosities η and ζ. The (total) energy density s=kBln[(e/n+ǫv0n)d/2(1/v0n 1)]+const. (4) eT =eˆ+ρv2/2 including the kinetic part is governedby − 1 ∂ eT = [eTv+(↔Π ↔σ) v]+ (λ T), (9) L = 202ℓ. Time is measured in units of t0 = mℓ2/ηv0. ∂t −∇· − · ∇· ∇ Here we assume constant viscosities with η = ζ. The strength of dissipation is represented by the normalized where λ is the thermal conductivity. We construct the stress tensor Π such that the entropy production rate viscosityη∗ =v0η/ℓ√mǫ. Wesetη∗ =0.11/2. Onallthe ij boundary walls we assume [13,14] within the fluid assumes the usual form [8], ∂ 1 n=nw+ℓwν n, (14) = dr [ λ T + σ v ]. (10) ·∇ ij i j ∂tS Z T ∇· ∇ ∇ Xij where ν is the outward normal unit vector at the sur- If there is no heat flow from outside, the above quan- face. We set nw = 5nc/2 and ℓw = 2√2ℓ, for which the wall is wetted by liquid in equilibrium in our simu- tity becomes non-negative-definite. To obtain Eq.(10) lation. We prepared an initial equilibrium state with a we need to require the general relation, circular gas bubble with radius 50 placed at the center 1↔ δ δ of the cell. The temperature was at T =0.875Tc, where Π = n eˆ . (11) ∇·(cid:18)T (cid:19) − ∇(cid:18)δnS(cid:19) − ∇(cid:18)δeˆS(cid:19) the density is 0.37nc in gas and 1.74nc in liquid. The eˆ n bottom boundary (y = 0) was then increased by a con- Using the van der Waals pressure p=nkBT/(1 v0n) stant ∆T, while the top boundary (y = L) was held at ǫv0n2 we then find − − the initialtemperatureTt =0.875Tc. (Hereweuse ”bottom”and”top”thoughnogravityisassumed.) Theside Πij =(p+p1)δij M in jn, (12) boundaries (x = 0 and L) are insulating (ν T = 0). − ∇ ∇ ·∇ Furthermore,thedensity-dependenceofthethermalcon- where p1 =M(n 2n+ n2/2)+Tn n (M/T). ∇ |∇ | ∇ ·∇ ductivity is crucial away from criticality [15], so we set As an example, we consider a bubble with radius R inheatflow. Thermocapillarybubble migrationhaslong λof=ord(ve0rη1/k(Binmu)nnitwsiothf ℓt2h/et0th=ervm0ηa/lmdi)ff.uTsihveintyλDiTn lbiqeuinidg been observed especially in space [9–11]. Theoretically, is larger than that in gas by 5 times. however,thebubblemotionhasbeenexaminedintheab- (i) First, for very small ∆T, the bubble position is sence of first-orderphase transition at the interface such shifted towardthe bottom without shape changes. Fig.1 as in the case of air bubbles in silicone oil [9]. In such shows such a steady state at ∆T = 0.00675Tc taken at casesthetemeraturedependenceofthesurfacetensionγ t = 3.4 103 after application of heat flux. Here the causes a bubble velocity of order vY = Rγ1(dT/dz)∞/η × migrationis stopped by the repulsionfromthe walls due in a given temperature gradient (dT/dz)∞ (Marangoni to the wetting condition in Eq.(14). The velocity field effect), where γ1 = −dγ/dT > 0 and η is the viscosity v is almost parallel to the interface, but there is a small in liquid. In one-comonent fluids the effect is drastically velocity component through the interface. Thus, at the alteredduetolatentheatgenerationorabsorptionatthe interface, first-order phase transition takes place and la- interface [12]. We point out that a convectivevelocity of ′ tent heat is generated or absorbed. As demonstrated in order vO =λ(dT/dz)∞/ρT∆s inside the bubble is suffi- ′ Fig.2, latent heat transport is so efficient that the tem- cient to carry the applied heat flux, whereρ is the mass perature gradient almost vanishes within the bubble. density in gas, ∆s is the entropy difference between gas (ii) Second, for the case ∆T = 0.054Tc, Fig.3 shows and liquid per unit mass, and λ is the thermal conduc- droplet migrationtowardthe bottom. Fig.4 displays the tivity in liquid. Note that the ratio velocity field in the steady state taken at t = 8 104. × M1 =vY/vO =Rγ1ρ′T∆s/ηλ (13) The velocity component through the interface is more increased than in Fig.1, again resulting in a flat temper- is a very large number of order R/a with a being a mi- ature inside the droplet as in Fig.5. As can be seen in croscopiclength. With phase changethe bubble velocity Figs.4 and 5, a thin liquid layer is sandwiched between is of order vO and is very slow [12]. Furthermore, the the bottom boundary and the droplet. The thickness of temperature inside the bubble shouldbe nearly homoge- this layer is so thin that we can define an apparent con- neousduetothelatentheatconvection. Thisisanatural tact angle θeff, which is a decreasing function of ∆T. In consequence in one-component fluids without contami- accord with these results, Garrabos et al. observed that nation[1,12],becausethe pressureshouldbecomenearly gas spreads on a heated wall initially wetted by liquid homogeneous outside the bubble and the fluid near the and exhibits an apparent contact angle [1]. interface should be close to the coexistence curve in the (iii) Third, we increased ∆T to 0.10125Tc, where the p-T phase diagram. bottom temperature is 0.97Tc. As shown in Fig.6, phase In our simulation we set K = 0 and then the inter- separationoccurredto produce a gasdomainatthe bot- face thickness is ℓ = (Cv0/2kB)1/2 far from the critical tom, into which the suspended droplet was finally ab- −1 point (since (µ µcx)/kBT is of order 1 for n v0 ). sorbed. In the steady state in Fig.7, the bottom is cov- − ∼ Space is measured in units of ℓ and the system length is ered by a gas layer or is completely wetted by gas [16], but the interface with the bulk liquid is still curved and ity in Physics) from the Ministry of Education, Culture, a small velocity field is induced. As can be seen in Fig.8 Sports, Science and Technology of Japan. taken at t = 3 104, the temperature gradient within × the gas layer becomes steeper than that in the liquid. It iswell-knownthatagasfilmona heatedwallsuppresses heat transfer to the bulk liquid and, as it thickens, film boiling is induced in gravity [3]. In the steady states the velocity is of order vO ∼ [1] Y. Garrabos, C. Lecoutre-Chabot, J. Hegseth, V. S. D ∆T/TL in terms of D = λ/ρC from C ∆s. AsTa result, the Reynolds Tnumber(∼ pPr−1R∆Tp/∼LT) is N64ik,o0l5a1y6ev02, D(2.0B01e)y.sens, and J.-P. Delville Phys. Rev. E verysmall,wherethePrandtlnumberPr =v0η/mDT is [2] L. M. Pismen and Y. Pomeau, Phys. Rev. E 62, 2480 here taken to be of order 1. In fact, in Figs.1, 4, and 7, (2000). the variance v2 is0.12,1.9,and0.63,respectively,in [3] V. S. Nikolayev and D. A. Beysens, Europhys. Lett, 47, units of 10−3pℓ/ht0,iwhere ℓ/t0 =DT/Prℓ. In our simula- 345 (1999). tion the viscous dissipation is very strong even in tran- [4] A. Onuki,Physica A 314, 419 (2002). sient states, while in the theory in Ref. [11] the inertia [5] A. Onuki, Phase Transition Dynamics (Cambridge Uni- effect (vapor recoil) is important in bubble motion on a versity Press, Cambridge, 2002). [6] J.S. Rowlinson, J. Stat. Phys. 20, 197 (1979). heatedwall. Thisdifferencestems fromthe factthatour droplet size is very small if ℓ is microscopic. [7] ForthevanderWaalsmodelγ ∝M1/2(1−T/Tc)3/2 well holds for T >0.8Tc within a few percents. Inthediffuse-interfacemethodsforcrystalgrowth[17] [8] L.D. Landau and E.M. Lifshitz, Fluid Mechanics (Perg- and thermo-capillary flow [18,19], the interface width ℓ amon, 1959). may be conveniently treated as a free parameter. This [9] N. O. Young, J. S. Goldstein, and M. J. Block, J. Fluid is because the phase field equations yield the moving- Mech. 6, 350 (1959). interface problems with appropriate interface boundary [10] B. Braun, Ch.Ikier,H.Klein andD.Woermann, Chem. conditions as long as ℓ R(= domain size). Also in Phys. Lett. 233, 565 (1995). ≪ fluids we may take ℓ as a free parameter in the range [11] D. Beysens, Y. Garrabos, V. S. Nikolayev, C. Lecoutre- ℓ R as far as we treat the steady state profiles (where Chabot, J.-P. Delville, and J. Hegseth, Europhys. Lett., ≪ 59, 245 (2002). the Marangoni effect is suppressed and the Stokes ap- [12] A.Onuki,tobepublished.Thetemperaturegradientin- proximation is valid). However, to describe rapid fluid ∗ sideadropletissmallerthatthatfarfromitby1/M1M2, motions, the normalized viscosity η (∝ η/ℓ) should be where M2 =Rρ′∆s/γ1 is also of order R/a@≫1. made very small if ℓ is enlarged. [13] P.G. deGennes, Rev. Mod. Phys. 57, 827 (1985). In summary, we have developed a scheme to study [14] S. Puri and K. Binder, Z. Phys.B, 86, 263 (1992). phase transition dynamics with inhomogeneous temper- [15] R. B. Bird, W. E. Stewart, and E. N. Lightfoot, Trans- ature under mass and heat fluxes. It can serve as a port Phenomena (Wiley, NewYork, 2002), p.272. phase field model in fluids, though we should examine [16] At large heat flux the density steeply changes within a its sharp-interface limit [19]. For one-component fluids few lattice pointsat theheated wall. Heretheboundary wementionwetting dynamicswithevaporationandcon- condition Eq.(16) is questionable. densation [2,13], as exemplified here, and boiling or con- [17] G. Caginalp, Phys. Rev. A 39, 5887 (1989); R. Kobayashi, Physica D 63, 410 (1993); A. Karma and vection in gravity [4,5]. These effects can in principle be W. J. Rappel,Phys. Rev.E 57, 4323 (1998). studiedinthe scheme ofthermocapillaryhydrodynamics [18] D. Jasnow and J. Vinãls, Phys. Fluids 8, 660 (1996); supplementedwithinterfaceconditions. However,hydro- R. Borcia and M. Bestehorn, Phys. Rev. E 67, 066307 dynamiccalculationsaccountingforlatentheattransport (2003). have rarely been performed. [19] D.M. Anderson, G.B. McFadden, and A.A. Wheeler, This work is supported by Grant for the 21st Cen- Annu.Rev.FluidMech.30,139(1998).Hereanattempt tury COE project (Center for Diversity and Universal- of taking thesharp-interface limit for fluidsis given. 3 FIG.1. Suspendedgas droplet in a steady state in heat flow with ∆T =0.00675Tc. The arrows indicate the velocity. ∆ T /T =0.00675 c T/T c 0.882 0.881 0.88 0.879 0.878 0.877 0.876 0 0.875 0.1 0.2 0.3 0.4 0.5 0 0.1 0.6 0.2 0.3 0.7 x 0.4 0.5 0.8 y 0.6 0.7 0.9 0.8 0.9 1 1 FIG.2. The temperature in the steady state in Fig.1. It is flat within thegas droplet due tolatent heat transport. 4 FIG. 3. Migration of a gas droplet for ∆T = 0.054Tc. It apparently wets the bottom partially in the steady state (left), where a thin liquid layeris still sandwiched. FIG.4. The velocity field around thegas droplet at t=8×104 in Fig.3. ∆ T /Tc= 0.054 T/Tc 0.93 0.92 0.91 0.9 0.89 0.88 0.87 0 0.1 0.2 0.3 0.4 0.5 0 0.1 0.2 0.3y 0.4 0.5 0.6 0.7 0.8 0.9 1 1 0.9 0.8 0.7 0.6x FIG.5. The temperature in the nearly steady state in Fig.4. 5 FIG.6. Time evolution for ∆T =0.10125Tc. Phase separation occurs at the bottom. The bottom is completely wetted by gas in the steady state. ∆T/T = 0.10125 c FIG.7. The velocity field around thegas region at t=3×104 in Fig.6. ∆ T /Tc= 0.10125 T/Tc 0.96 0.94 0.92 0.9 0.88 0 0.1 0.2 0.3 0.4 0.5 0 0.1 0.2 0y.3 0.4 0.5 0.6 0.7 0.8 0.9 1 1 0.9 0.8 0.7 0.6 FIG.8. The temperature in the nearly steady state in Fig.7. 6