Effects of Temperature Difference and Heat Loss on Oscillation Characteristics of Thermo-Solutocapillary Convection in Toluene/N-Hexane Mixed Solution

Abstract

During the crystal growth process using the floating zone method, the uneven distribution of impurities on the surface of the melt can trigger a coupling mechanism between solutocapillary convection driven by the concentration gradient and thermocapillary convection driven by the temperature gradient, resulting in the Marangoni convection at the free surface. When the temperature and concentration gradients reach certain values, the crystal surface and interior exhibit time-dependent, periodic oscillations, leading to the formation of micrometer-scale impurity stripes within the crystal. This study focuses on the effects of temperature difference and heat loss in a liquid bridge under microgravity on the structure and interface oscillation characteristics of thermo-solutocapillary convection, aiming to explore the coupling phenomenon of this oscillation and provide valuable information for crystal growth processes. An improved level set method is employed to accurately track every displacement of the interface, while the surface tension is addressed using the CSF model. In addition, the area compensation method is used to maintain simulation quality balance. A comprehensive analysis is performed on the oscillation characteristics of thermo-solutocapillary convection at the free surface, ranging from the temperature, concentration, deformation, and velocity distributions at the upper and middle heights of the liquid bridge. The results indicate that under small temperature differences (ΔT = 1 − 3), the transverse velocity at the upper end exhibits a single-periodic oscillation, while the longitudinal velocity presents a double-periodic oscillation. At the intermediate height, both the transverse and longitudinal velocities display a single-periodic oscillation. Under a large temperature difference (ΔT = 6), the oscillation of velocities at the upper end and the middle position become multi-periodic. In addition, heat loss has certain regular effects on the oscillatory flow of thermo-solutocapillary convection within a certain range. The velocity, amplitude, and frequency of the upper end and the middle position at the free surface decrease gradually, and the oscillation intensity also weakens with the increase in heat loss (Bi = 0.2 − 0.6). These new discoveries can provide a valuable reference for optimizing the crystal growth process, thereby enhancing the quality and performance of crystal materials.

1. Introduction

Marangoni convection is very common in many industrial applications, such as crystal growth [1,2,3], welding [4,5,6,7], and evaporation [8,9]. During crystal growth, the buoyancy convection is weakened under microgravity, and Marangoni convection becomes the dominant convection. The prepared crystals are not only large in size but also uniformly distributed [10]. The floating zone method is a crucial method for producing high-quality semiconductor materials. During the floating zone’s crystal growth, the semiconductor materials in a molten state are frequently doped with impurities. The non-uniform distribution of impurities on the free surface of the melt can generate the solute capillary convection due to the concentration gradient. The convection is coupled with thermocapillary convection driven by the temperature gradient, and they form Marangoni convection on the free surface eventually [11,12]. A liquid bridge model can be utilized to explore this phenomenon. The liquid bridge model is a model in which a fluid medium is located between two coaxial solid disks, and the upper and lower disks experience various temperature and concentration differences, respectively, thereby forming the thermo-solutocapillary convection on the surface of the liquid [13,14]. After the upper disk is heated and a certain amount of solute introduced, the surface tension based on temperature and concentration begins to change and forms Marangoni stress at the free interface, resulting in the fluid movement at the free interface. In situations where the disparity in temperature and concentration between the solid disks is small, a steady axisymmetric vortex structure will appear inside the liquid bridge. Simultaneously, the time-based periodic oscillation occurs internally when the temperature and concentration differences reach critical values. Chun et al. [15] observed this phenomenon through experiments for the first time.

During the floating zone technique for crystal growth, the primary factor responsible for generating micron-level impurity fringes within the crystal is oscillation [16]. Many scholars applied magnetic field [17,18], transverse vibration [19,20,21,22], and geometric size adjustment [23,24,25] to investigate the internal flow and oscillatory behavior of the thermo-solutocapillary convection aiming to suppress the oscillation. Takagi et al. [26] first studied the transition state of the liquid bridge from the static to oscillation at a low Prandtl number. The non-contact detection method was used in the experiment to premix tin oxide particles in the fluid. The results exhibit that the temperature begins to oscillate periodically and the corresponding oscillation frequency and amplitude are obtained when the Marangoni number increases to 43.3 in a high vacuum environment. Simultaneously, the frequency and amplitude of oscillation increase when the Marangoni number reaches 91. Wang et al. [27] investigated the flow structure of buoyant thermocapillary convection and the critical condition of oscillation under diverse height–diameter and volume ratios in ground experiments. Jayakrishnan and Tiwari [28] carried out a numerical research of oscillating thermocapillary flow in a 5cSt silicone oil-filled liquid enclosure. Under microgravity conditions, the finite volume method was applied to calculate in the liquid bridge with various volume ratios to determine the dynamic evolution of temperature disturbance during the beginning of oscillation instability. Zhou et al. [29,30] conducted a numerical study on the thermocapillary convection of nanofluids with medium and large Prandtl numbers in a rectangular cavity. The critical temperature difference of the unsteady flow is obtained via varying the thermal gradient, and the influence of the volume fraction of nanoparticles on the critical temperature gradient is analyzed. Finally, the corresponding oscillation laws are obtained. Dressler and Sivakumaran [31] verified through experiments that the stability of the liquid bridge is significantly affected by the heat loss of the bridge. Kamotani et al. [32] examined the impact of free surface heat loss on the oscillatory regime of thermocapillary-driven convection in a high Prandtl number liquid bridge via experimental and numerical analyses. They concluded that the heat loss effect is caused by the interaction between the dynamic surface morphology of the liquid bridge and the heat loss from the hot wall. Wang et al. [33] conducted experimental investigations on the oscillatory flow of thermocapillary convection in a liquid bridge and numerically analyzed the impact of heat loss on the initial stage of surface oscillation.

Thermo-solutocapillary convection is more complicated than pure thermocapillary convection. Since Bergman [34] first reported the thermo-solutocapillary convection, the study of thermo-solutocapillary convection has attracted more and more attention. Zhou and Huai [35] paid attention to the dynamic free surface of thermo-solutocapillary convection. The numerical results exhibit that the thermo-solutocapillary convection with various capillary ratios displays three distinct modes of free surface shape at low Marangoni values (small temperature difference). Hirata et al. [36] discussed the impact of the thermal and solute Marangoni numbers on Czochralski growth under microgravity conditions. Witkowski and Walker [37] investigated how the concentration distribution in SiGe alloy changes during growth, as a function of the growth rate. The concentration gradient no longer maintains a uniform distribution on the free interface when the growth rate (0.1 < V_g < 4) is considered, and the initiation of flow is primarily attributed to the concentration gradient present at the melting interface which creates a diffusive boundary layer. In this case, the instability of the growth zone comes from the flow inertia and the interface concentration gradient. The instability is not limited to strong concentration regions, but spreads to the whole melt. Zou et al. [38] studied the coupling of thermo-solutocapillary convection in the SixGe1-x system of a liquid bridge of a half floating zone by numerical simulation. The impact of the aspect ratio on the stability of capillary convection is also studied. Alhashash and Saleh [39] studied the coupling of solute and thermocapillary convection under buoyancy via numerical methods. The law of heat and mass transfer rate variation is obtained within a certain range of thermal Marangoni numbers (a certain temperature). Sarma and Mondal [40] studied Marangoni instability in binary mixtures heated from the free surface, taking into account the impact of thermocapillary convection, solutocapillary convection, and buoyancy. They analyzed the stability characteristics of the polymer located at the deformable interface and the rigid substrate with poor conductivity. Meanwhile, they concluded that the system has two different kinds of oscillatory instability besides the monotone perturbation. Fan and Liang [41] considered the Soret effect and studied the influence of the thermal Marangoni number on thermo-solutocapillary convection in a liquid bridge containing a solution of n-decane and n-hexane under microgravity by a numerical method. The results show that the thermo-solutocapillary convection wave occurs when the Ma_T number is large (large temperature difference), and the Ma_T number has an impact on the flow field structure of convection oscillation in the liquid bridge. Moreover, the fluctuation of thermo-solutocapillary convection increases with the increase in Ma_T number. Lyubimova and Skuridin [42] numerically investigated the feasibility of axial vibration on the evolution of convection in a liquid bridge under the neglect of gravity. It was observed that the convection varied with temperature when the solute Marangoni number remained constant. Furthermore, a time lag was found in the convection process, and simultaneous vortices were identified that were dominated by thermocapillary and solutocapillary effects, respectively. Zhuang and Zhu [43] employed numerical simulation techniques to examine the impact of buoyancy on the thermo-solutocapillary convection of non-Newtonian fluids. Their results indicated that as buoyancy increased, and the influence of surface tension on heat and mass transfer decreased. Furthermore, they observed an increase in both the average Nusselt number and Sherwood number with a rise in temperature, which was attributed to the effects of the Marangoni number and thermal Rayleigh number.

In summary, temperature difference in the thermo-solutocapillary flow refers to the presence of a temperature gradient on the liquid bridge surface, driving the flow of fluid near the surface and affecting the crystal growth process. Heat loss represents the energy loss that occurs on the liquid bridge surface, which may impact the temperature and concentration fields driving the flow, leading to changes in the oscillation characteristics of the flow. The oscillation characteristics of thermo-solutocapillary flow refer to the periodic or stochastic fluctuation of temperature, concentration, velocity, and free surface deformation within the liquid bridge. These oscillations are due to the interaction between temperature–concentration convection, gravity, surface tension, and buoyancy, among other factors. However, studies on the coupled thermo-solutocapillary convection with different temperature differences and heat losses within the float zone are still scarce, particularly considering dynamic interface deformation. Numerical simulations with accuracy and stability can be difficult to achieve under these circumstances. Therefore, the purpose of this study is to investigate these issues, taking into account the effects of the dynamic free surface deformation, and offer assistance in the crystal growth process using the float zone method. Section 2 of the study introduces the physical and mathematical models of thermo-solutocapillary convection and the level set method. Section 3 discusses the influence of different temperature differences on the convection structure and interface oscillation characteristics of thermo-solutocapillary flow. Subsequently, we discuss the impact of heat loss on the convection structure and oscillation characteristics of thermo-solutocapillary flow. Finally, the effects of different variables on the flow and oscillation at the liquid bridge interface are summarized.

2. Models and Methods

2.1. Physical Model

Figure 1 illustrates the physical model of the liquid bridge that initially forms between two horizontal plates. The liquid bridge has an initial length of R and a height of H and is surrounded by a gas medium in a rectangular container with a height of H and a length of 2R. The liquid bridge comprises a binary mixture of toluene and n-hexane with a concentration ratio of 0.24 to 0.76. The initial state of the gas–liquid interface is vertical. The upper and lower plates are subjected to various temperatures and concentrations that satisfy T₁ > T₀ and C₁ > C₀, that is, the temperature difference between the upper and lower plates can be expressed as ΔT = T₁ − T₀, while the concentration difference is represented by ΔC = C₁ − C₀.

A succession of assumptions is introduced to simplify the research. Firstly, it is assumed that the upper and lower plates are fully insulated and impermeable. Additionally, the slip-free condition is applied to all walls of the computational domain. The impact of viscous dissipation is considered negligible, and it is assumed that the flow is laminar. Moreover, gravitational effects are ignored. Finally, the binary mixture fluid under investigation is an incompressible Newtonian fluid. It should be noted that the physical parameters of the fluid remain unchanged, except for the fluid surface tension. The fluid’s surface tension can be expressed as a linear function with respect to both temperature and concentration, which is interpreted as σ = σ₀ − σ_T (T − T₀) − σ_C (C − C₀), where σ₀ represents the initial surface tension, T and C represent coefficients of the temperature and concentration of the fluid, respectively, and σ_T and σ_C represent the coefficient of surface tension varying with temperature and concentration, respectively. Furthermore, σ_T = −∂σ/∂T = constant and σ_C = −∂σ/∂C = constant.

2.2. Governing Equations and Boundary Conditions

This study mainly investigates the thermo-solutocapillary convection under microgravity conditions. To ensure the precision of the final outcomes, the surface deformation is considered. This study uses the level set method to capture the interface displacement in the application of a gas–liquid two-phase system.

To facilitate computational efficiency, it is imperative that the governing equations be transformed into a dimensionless form. The characteristic density and viscosity are ρ_l and μ_l. The characteristic length is L = H = 0.5R, and the characteristic velocity is U = (σ_T·ΔT)/μ_l. The subscripts g and l represent gas and liquid, respectively. The ensuing expression encapsulates the dimensionless parameters designated:

$$\overline{x}=\frac{x}{L}, \overline{z}=\frac{z}{L}, \overline{u}=\frac{u}{U}, \overline{v}=\frac{v}{U}, \overline{\rho }=\frac{\rho }{{\rho }_l}, \overline{\mu }=\frac{\mu }{{\mu }_l}, \overline{p}=\frac{p}{{\rho }_l{U}^2}, \overline{t}=\frac{t}{L/U}$$

(1)

where the physical quantity with the “-” symbol as the superscript represents a dimensionless quantity. $\overline{x}$ and $\overline{z}$ are utilized to represent dimensionless distances in the horizontal and vertical directions, respectively. $\overline{u}$ and $\overline{v}$ refer to dimensionless transverse and longitudinal velocities, respectively. $\overline{\rho }$, $\overline{\mu }$, $\overline{p}$, and $\overline{t}$ refer to dimensionless density, viscosity, pressure, and time, respectively.

The governing equations consist of the following dimensionless mass conservation, Navier–Stokes, energy conservation and concentration diffusion equations, and interface capture function [35], as listed in Formulas (2)–(6). Among them, the energy conservation and concentration diffusion equations are derived based on the continuum assumption and the conservation of mass, momentum, and energy [44], using the equation of state and Fick’s law of diffusion, respectively.

$$\nabla \cdot \mathbf{V}=0$$

(2)

$$\begin{array}{ll}\overline{\rho }\left[\frac{\partial \mathbf{V}}{\partial \overline{t}}+\nabla \left(\mathbf{V}\cdot \mathbf{V}\right)\right]& =-\nabla \overline{P}+\nabla \left(\frac{2\overline{\mu } {\mathrm{D}}^{\prime }}{Re}\right)+\frac{\left(C{a}_T{\nabla }_S\Theta +C{a}_C{\nabla }_S{C}^{\prime }\right)}{We}\\ & +\frac{\left[\left(1-C{a}_T\Theta -C{a}_C{C}^{\prime }\right)\cdot \kappa \left(\phi \right)\delta \left(\phi \right)\nabla \phi \right]}{We\cdot \overline{\rho }}\end{array}$$

(3)

$$\frac{\partial \Theta }{\partial \overline{t}}+\mathbf{V}\cdot \nabla \Theta =\frac{{\nabla }^2\Theta }{Ma}$$

(4)

$$\frac{\partial C^{\prime }}{\partial \overline{t}}+\mathbf{V}\cdot \nabla C^{\prime }=\frac{{\nabla }^2{C}^{\prime }}{\left(Ma\cdot Le\right)}$$

(5)

$$\frac{\partial \phi }{\partial \overline{t}}+\mathbf{V}\cdot \nabla \phi =0$$

(6)

where $\mathbf{V}=(\overline{u},\overline{v} )$ represents the dimensionless fluid velocity, D′ represents the viscous stress tensor, κ represents the interface curvature, δ represents the Dirac delta function, φ represents the level set function, and ∇_s = (I − nn)·∇ represents the free interface gradient operator. Θ and C′ refer to the dimensionless temperature and concentration, respectively, taking Θ = (T − T₀)/(T₁ − T₀) and C′ = (C − C₀)/(C₁ − C₀). Moreover, the dimensionless quantity is expressed as the Reynolds number Re = (ρ_lUH)/μ_l, Weber number We = (ρ_lU²H)/σ, Prandtl number Pr = μ_l/(ρ_lα_l), Marangoni number Ma = Re·Pr, thermal capillary number Ca_T = (σ_T·ΔT)/σ₀, solutal capillary number Ca_C = (σ_C·ΔC)/σ₀, and Lewis number Le = α_l/D_l, where α_l and D_l refer to the thermal and solutal diffusivity, respectively.

The gas phase boundary conditions are divided into adiabatic and heat loss. The condition of zero gradient is imposed on all walls. The boundary and initial conditions can be expressed as:

Initial time t = 0:

$$\Theta =0, C^{\prime }=0, \overline{p}=0, \mathbf{V} =0, \phi \left(1,Z\right)=0$$

(7)

Upper plate (Z = H):

Θ = 1, C′ = 1, V = 0, ∂φ/∂Z = 0

(8)

Lower plate (Z = 0):

Θ = 0, C′ = 0, V = 0, ∂φ/∂Z = 0

(9)

Adiabatic: Outer edge of gas phase (X = 0 and X = 4),

V = 0, ∂Θ/∂X = 0, ∂ C′/∂X = 0, ∂φ/∂X = 0

(10)

Heat loss: Outer edge of gas phase (X = 0 and X = 4),

V = 0, ∂Θ/∂X = −Bi·Θ, ∂ C′/∂X = 0, ∂φ/∂X = 0

(11)

where Bi = h·L/λ_l, h is the convective heat transfer coefficient.

2.3. Level Set Method

The level set method is a numerical technique used to represent and track moving interfaces between two fluids. In this method, the material interface is represented as a zero level isosurface of the function φ (x, t), which is a level set function. The level set function has a positive value outside the liquid bridge, a zero-value interface between gas and liquid and a negative value inside the liquid bridge. By solving a convection equation using the level set function, the motion of the interface can be predicted and tracked over time. This method was first introduced by Osher and Sethian [45] and has since become widely used in various computational applications related to fluid dynamics and material sciences. With respect to the level set function φ (x, t), the movement of the interface is forecasted by solving the subsequent convection equation:

$${\phi }_t+\mathbf{V}\cdot \nabla \phi =0$$

(12)

To ensure the level set function φ retains its status as a distance function, it is necessary to perform a process known as level set function reinitialization. This process involves updating the level set function to conform to the expected properties of a distance function [46]:

$$\{\begin{array}{c}{\phi }_t=\mathit{sign}\left({\phi }_0\right)\left(1-\left|\nabla \phi \right|\right)\hfill \\ \phi \left(X,0\right)={\phi }_0\left(X\right)\begin{array}{cc}& \end{array}\hfill \end{array}$$

(13)

where $\left|\nabla \phi \right|=\sqrt{{\phi }_x^2+{\phi }_{\mathrm{z}}^2}$. To simplify the solution, the sign(φ₀) needs to be smoothed as [14]:

$$sign({\phi }_0)=\phi /\sqrt{{\phi }^2+{\epsilon }^2}$$

(14)

Despite the completion of the re-initialization process of the level set function, the area conservation cannot be reached yet due to the numerical dissipation. Therefore, it is recommended to solve the equation for the compensation of the area until it reaches a stable state [47]:

$$\frac{\partial \phi }{\partial \overline{t}}+\left[1-A\left(t\right)/A_0\right]f\left(c\right)\left|\nabla \phi \right|=0$$

(15)

where A(t) represents the liquid bridge area at time t as determined by the level set function, A₀ represents the initial area of the liquid bridge. The function f(c) is a representation of the area constraint and can be regarded as a function of the local curvature. It is subject to variation with respect to the parameters h and n, which is expressed as:

$$f(c)<1.0+h^{\prime }{\gamma }^n$$

(16)

The values h′ = 0 and n = 0 are denoted, as they result in an accelerated convergence of the numerical procedure. The criteria for convergence can be formulated as follows:

$$\left|A_0-A\left(\overline{t}\right)\right|<1.0\times 10^{-4}$$

(17)

By employing this method, the total mass is conserved with a high level of accuracy, not deviating beyond 0.01%. To mitigate against numerical instability on the interface, particularly in the event of a significant density ratio existing within the two-phase system, it becomes necessary to smooth the physical properties in close proximity to the interface. The interface is assumed to comprise of a thin layer that possesses a thickness of 2α and maintains continuity without abrupt alterations concerning physical properties. To aid in this process, both the Heaviside and Dirac functions are introduced, and their definitions are listed as follows [48]:

$$H_{\alpha }\left(\phi \right)=\{\begin{array}{c}0, if \phi <-\alpha \hfill \\ \left[1+\phi /\alpha +\sin \left(\pi \phi /\alpha \right)\right]/2, if \left|\phi \right|\le \alpha \hfill \\ 1, if \phi >\alpha \hfill \end{array}$$

(18)

$${\delta }_{\alpha }\left(\phi \right)=\{\begin{array}{l}\left(1+\cos \left(\pi \phi /\alpha \right)\right)/2\alpha , if \left|\phi \right|<\alpha \\ 0, if \left|\phi \right|\ge \alpha \end{array}$$

(19)

The density ρ and dynamic viscosity μ can be exhibited as follows:

$$\{\begin{array}{c}\rho \left(\phi \right)={\rho }_l+\left({\rho }_g-{\rho }_l\right)H_{\alpha }\left(\phi \right)\\ \mu \left(\phi \right)={\mu }_l+\left({\mu }_g-{\mu }_l\right)H_{\alpha }\left(\phi \right)\end{array}$$

(20)

A uniform staggered grid is employed to discretize the liquid bridge and the surrounding air region, using a central difference scheme for all derivatives except the convection term. The QUICK method is used for convection terms in the Navier–Stokes equation to achieve third-order accuracy, while the central difference scheme is adopted for other terms with first-order accuracy. The level set function equation is evaluated using the QUICK method for third-order accuracy. The Poisson equation for pressure is solved using the successive over relaxation (SOR) method, while the continuum surface force (CSF) model is utilized to handle surface tension at the interface. Additionally, the ENO scheme is implemented to solve the convection term of the level set function equation. Overall, the methods employed in this simulation exhibit high accuracy, and the adopted grid discretization and difference scheme are capable of accurately capturing the characteristics of the liquid bridge and surrounding air.

The liquid bridge working medium is a toluene/n-hexane mixed solution, and physical properties [49,50] are exhibited in

Table 1. Physical properties of toluene/n-hexane mixture at 25 °C (0.24/0.76).

Property	Symbol	Value
Density	ρl	699 kg/m3
Diffusion coefficient	αι	1.88 × 10−7 m2/s
Dynamic viscosity	${\mu }_{\iota }$	3.37 × 10−4 kg/(m·s)
Surface tension	σ0	2.1 × 10−2 N/m
Temperature surface tension coefficient	${\sigma }_T$	9.43 × 10−5 N/(m·K)
Concentration surface tension coefficient	${\sigma }_C$	−8.62 × 10−3 N/m
Prandtl number	Pr	5.54
Lewis number	Le	25.8

.

This study applies a variety of different grids for detection calculations to ensure the grid independence. The report of grid refinement is exhibited in

Table 2. Dimensionless velocity, temperature, and concentration distributions at the intermediate height surface point of the right liquid bridge (ΔT = 1, ΔC = 1, Pr = 5.54, Le = 25.8, R = 5 mm, Ar = 1.0).

Grids	Longitudinal Velocity	Transverse Velocity	Temperature	Concentration
41 × 21	0.005382	−0.000136	0.265457	0.316837
81 × 41	0.005194	−0.000129	0.265165	0.313024
101 × 51	0.005234	−0.000131	0.265960	0.314206
121 × 61	0.005162	−0.000127	0.264634	0.312122

. To validate the code and refine the mesh, we selected the velocity, temperature, and concentration distributions at the intermediate height interface of the liquid bridge as the objects of analysis. Based on the error calculation method utilized, the values for error in longitudinal velocity comparison are 3.6%, 0.7%, and 0.6% for the second group and other groups, respectively. Correspondingly, the error values for the transverse velocity comparison are 5.4%, 1.6%, and 1.6%. Moreover, in the temperature comparisons, the error values are 1.1%, 0.3%, and 0.2%, while in the concentration comparisons, the respective values are 1.2%, 0.4%, and 0.3%. Following the analysis of the results, we concluded that a square grid measuring 81 × 41 is sufficient to capture the interface, temperature, concentration, and flow fields with a high degree of accuracy.

In order to verify the present model and code, we compared present results with Zhou’s results [35]. Other calculation parameters are R_ω = −10, Le = 100, Ar = 1.2, and Pr = 1. Figure 2 shows the surface longitudinal velocity and surface deformation in different Ma numbers, respectively. The present results are consistent with Zhou’s results, and these comparisons illustrate the correctness of the present model and code.

3. Results and Discussions

This paper mainly studies the influence of temperature difference and heat loss on the flow structure and surface oscillation of the coupled thermo-solutocapillary convection. The calculation conditions are exhibited in

Table 3. Calculation conditions (corresponding to Section 3.1).

ΔT(-)	ΔC(-)	CaT	CaC	Ma
1	1	0.00045	−0.021	372
2	1	0.0009	−0.021	744
3	1	0.00135	−0.021	1116
6	1	0.0027	−0.021	2232

and

Table 4. Calculation conditions (corresponding to Section 3.2).

ΔT(-)	ΔC(-)	CaT	CaC	Ma	Bi
1	1	0.00045	−0.021	372	0.2
1	1	0.00045	−0.021	372	0.4
1	1	0.00045	−0.021	372	0.6

, respectively. Two fixed monitoring points “a” and “b” are installed at the upper end (R, 0.9H) and the intermediate position (R, 0.5H) of the right interface to monitor the oscillation characteristics of the free interface.

3.1. Effect of Temperature Difference on Flow Structure and Interfacial Oscillation Characteristics of Thermo-Solutocapillary

Figure 3 shows the velocity vector under different temperature differences at the same time. It can be found in Figure 3 that the adjustment of temperature difference has little effect on the velocity vector in the vicinity of the central axis of the liquid bridge (x = 1.75 to 2.25), but the velocity vector in other areas has changed greatly due to the change in temperature difference. It can be clearly examined from Figure 3a–d that the flow velocity of the fluid near the interface on both sides of the liquid bridge escalates with an augmentation in temperature difference. The fluid flow near the interface is relatively stable when the temperature difference is 1 to 3. However, when the temperature difference increases to 6, the fluid flow near the interface is no longer stable and fluctuates intensively, resulting in the interface fluid flow becoming disordered, thereby promoting the absorption of liquid bridge surface materials. Because the interface vicinity is the principal surface recirculation area, the increase in the flow rate indicates that the conversion speed between volume reflow and surface reflow is accelerated. By comparing the velocity vector diagrams under different temperature differences, it indicates that the velocity of fluid in the vicinity of the interface increases with the increase in temperature difference. When the temperature difference is aggravated to a particular extent, the interface flow becomes disordered. The conversion speed of surface reflow and volume reflow also increases as the temperature difference increases.

Figure 4a exhibits the schematic diagrammatic sketch of the working area of the liquid bridge. The free surface on the right side of the liquid bridge is selected as the research object. Figure 4b exhibits the temperature change of the right free surface, and Figure 4c shows the concentration change of the right free surface. It can be observed from Figure 4b that the surface temperature is distributed as a convex function under different temperature differences. Additionally, a notable disparity in the temperature gradient is noticed between the upper and lower ends, with the former exhibiting a significantly larger variation. As the temperature difference between the upper and lower ends increases, the gradient difference between these two ends diminishes. This suggests that the temperature distribution at the interface becomes more uniform as the temperature difference increases. Concerning the concentration change of the free interface, it can be found from Figure 4c that the hot end (z = 0.85 to 1.0) and cold end (z = 0 to 0.15) of the liquid bridge present a large concentration gradient. The concentration increased sharply from 0.3 to 1 near the hot end and decreased sharply from 0.2 to 0 near the cold end. Because the hot and cold ends near the interface are the conversion areas of volume reflow and surface reflow, a large concentration gradient is formed in this area. Simultaneously, it can be found that in the range of height z = 0.15 to 0.85 the concentration shows a decreasing trend as the temperature difference increases. The temperature change presents more uniform than the concentration change.

Figure 5 exhibits the transverse and longitudinal velocities of the free surface on the right side, as affected by various temperature differences. It can be found from Figure 5a that the transverse velocity generally presents a transverse S-type distribution. The peak of transverse velocity under different temperature difference conditions is close to the hot end (z = 0.7 to 1.0) and the cold end (z = 0 to 0.3) of the system, and the velocity varies sharply. The velocity at the intermediate position (z = 0.3 to 0.7) does not vary significantly with the temperature difference. The peak values of transverse velocity appear near z = 0.15 and z = 0.85 near the cold and hot ends, respectively. The peak value rises with an elevation in the temperature gradient, while the velocity value at the middle height is relatively stable. The velocity value fluctuates greatly as the temperature difference increases. The conversion process between surface flow and volume reflux of the system is accelerated. The acceleration of surface flow will also lead to a big pressure gradient around the interface. The pressure gradient near the interface increases the transverse velocity fluctuation. As listed in Figure 5b, the longitudinal velocity presents a parabolic shape, and there is no complex fluctuation phenomenon. The longitudinal velocities are all positive, indicating a unidirectional fluid flow from the lower to the upper plate. Due to the high concentration of the hot end and the absorption of the interface solute, the interface has a higher concentration gradient. At this time, the solutocapillary force exceeds the thermocapillary force and dominates the orientation of the interfacial fluid motion.

Figure 6 exhibits the deformation trend of the free surface on the right side, as affected by various temperature differences. On the whole, the interface deformation is S-type. The interface deformation is convex from z = 0 to 0.5 and concave at z = 0.5 to 1. The degree of convexity and concaveness increases as the temperature difference increases. Because the near-hot end temperature of the system exceeds the near-cold end temperature, as well as the gas side pressure at the upper end surpassing the lower end’s gas side pressure, the overall surface deformation is concave and convex from top to bottom, respectively. Moreover, the conversion process between surface flow and volume reflow of the system will accelerate with an increase in the temperature differential. The internal disturbance of the liquid bridge will intensify, resulting in a rise in the convexity and concavity of the interface with the augmentation of the temperature differential.

Figure 7a exhibits the oscillation distribution of the transverse velocity at the upper monitoring point “a” for various temperature differences. The location of the monitoring point is listed in Figure 1, and the coordinate is (R, 0.9H). Figure 7b shows the power spectral density (PSD) of transverse velocity. As shown in Figure 7a, the transverse velocity oscillation at the upper monitoring point is approximately a small amplitude sinusoidal oscillation when the temperature differential is small in the yellow portion. The mean values of transverse velocity are −0.0025, −0.00276, and −0.00294, respectively. The amplitudes are 2.5 × 10⁻⁵, 3.4 × 10⁻⁵, and 4.0 × 10⁻⁵, respectively. The amplitude magnifies as the temperature increases. Furthermore, the magnitude of the velocity oscillation escalates from 10⁻⁵ to 10⁻³ when the temperature difference rises to 6. The oscillation is converted into a large amplitude pulsation oscillation. The oscillation amplitude is large and irregular. As listed in Figure 7b, the oscillation power density (in the yellow portion) at small temperature differences slowly increases as the temperature difference rises, ranging from 7.5 × 10⁻⁹ to 3.9 × 10⁻⁸. For sinusoidal oscillations under small temperature difference conditions, there is only one fundamental frequency in the corresponding PSD diagram. The frequencies are f₁ = 0.16, f₂ = 0.17, and f₃ = 0.13, respectively. Therefore, the oscillation under small temperature difference is a single-periodic oscillation. As for oscillations under a large temperature difference, there are multiple dominant frequencies, including f₄ = 0.014, f₅ = 0.035, f₆ = 0.1, f₇ = 0.11, and f₈ = 0.21. In the main frequency, the power density of the oscillation is between 1.2 × 10⁻⁵ and 2.0 × 10⁻⁵. These five frequencies are independent to each other, so the pulsating oscillation under large temperature differences is a multi-periodic oscillation.

Figure 8a exhibits the oscillation distribution of the longitudinal velocity of the upper monitoring point “a” under different temperature difference conditions, and Figure 8b is the PSD of longitudinal velocity. Figure 8a shows that the longitudinal velocity exhibits a small amplitude pulsation oscillation characteristic under the small temperature difference condition, and the oscillation exhibits a regular parabolic shape. The time interval for each parabolic shape is approximately 40 when the temperature difference is 2. The time interval for each parabolic shape is approximately 32 when the temperature difference is 3. The fact that the longitudinal velocities at different temperature differences are all positive suggests the presence of a unidirectional flow of fluid from the lower plate to the upper plate. The mean values of longitudinal velocity are 0.00131, 0.00164, and 0.00191, respectively. The amplitudes are 1.3 × 10⁻⁵, 1.5 × 10⁻⁵, and 1.6 × 10⁻⁵, respectively. The amplitude increases as the temperature difference rises. The average speed is 0.00224 when the temperature difference rises to 6. The amplitude increases exponentially, and the oscillation becomes a large amplitude pulsating oscillation. Because the upper monitoring point is in the active zone, the surface reflux intensifies the flow as the temperature difference rises, leading to an increase in the amplitude of longitudinal velocity. The oscillation evolves into a large amplitude pulsating oscillation when the surface flow accelerates to a certain extent. It can be found from Figure 8b that the spectrum of the small temperature difference has two main frequencies. The two corresponding fluctuations in the spectrum are independent of each other. Smaller fluctuations correspond to low-frequency oscillations, while larger fluctuations correspond to high-frequency oscillations. The small-amplitude pulsation oscillation at this time is the double-periodic oscillation. The oscillation energy also has little difference under different temperature differences, with power spectral densities ranging from 2.6 × 10⁻¹⁰ to 2.7 × 10⁻¹⁰, from 1.3 × 10⁻⁹ to 2.2 × 10⁻⁹, and from 9.7 × 10⁻¹⁰ to 2.0 × 10⁻⁹, respectively. Multiple main frequencies appear in the spectrogram under large temperature differences, including f₇ = 0.033, f₈ = 0.04, f₉ = 0.19, and f₁₀ = 0.2. In the main frequency, the power density of the oscillation is between 2.1 × 10⁻⁶ and 3.0 × 10⁻⁶. At present, the oscillation is characterized by the occurrence of multiple periodicities. Compared with the transverse velocity in Figure 7, the oscillations in the longitudinal velocity display a double-periodic oscillation under small temperature differences, which is different from the single-periodic oscillation of the transverse velocity. However, they are all multi-periodic oscillation under large temperature differences. Meanwhile, the power of the longitudinal velocity is smaller than the power of the transverse velocity, indicating that the velocity in the active region of the hot corner changes relatively strongly in the transverse direction.

Figure 9a shows the transverse velocity oscillation distribution of the monitoring point “b” under various temperature differences. The coordinate of the monitoring point is (R, 0.5H) which is in the middle position of the right interface. Figure 9b is the PSD of transverse velocity. As listed in Figure 9a, under conditions of minimal temperature variation, the average values of velocity are −0.000129, −0.000176, and −0.000238, respectively. The amplitudes are 2.5 × 10⁻⁶, 1.2 × 10⁻⁶, and 1.8 × 10⁻⁶, respectively. The amplitude also magnifies as the temperature difference rises. The oscillation curve gradually transforms to a regular large-amplitude oscillation from a small-amplitude pulsating oscillation as the temperature difference increases. This time interval of regular oscillation is approximately 64. The oscillation enters the large-amplitude pulsation phase when the temperature difference increases to 6. Because the intermediate position of the interface is closer to the vortex core, the transverse velocity oscillation is affected by the micro-scale lateral migration of the vortex core as the temperature difference increases. The lateral migration scale increases result in an increase in the amplitude of transverse velocity. Upon reaching a particular threshold of temperature difference, the vortex core will be severely reciprocating. The transverse velocity thus enters the large-amplitude pulsation oscillation phase. The oscillation of the middle monitoring point exhibits more pulsating oscillation characteristics than the transverse velocity oscillation of the upper monitoring point. Figure 9b shows that there is only one main frequency in the spectrum when there is a small temperature difference. The frequencies are f₁ = 0.242, f₂ = 0.244, and f₃ = 0.24, respectively. The power spectral density is between 2.8 × 10⁻¹¹ and 2.8 × 10⁻⁹ and slowly increases as the temperature difference increases. Therefore, both small-amplitude pulsation oscillations and regular oscillations under small temperature differences are single-periodic oscillations. Multiple main frequencies appear in the spectrum diagram when the temperature difference is 6. These frequencies are independent of each other, so the large-amplitude pulsation oscillation under a large temperature difference is the multi-periodic oscillation.

Figure 10a exhibits the oscillation of the longitudinal velocity of the monitoring point “b” under different temperature difference conditions, and Figure 10b is the PSD of longitudinal velocity. As listed in Figure 10a, the magnitude of oscillation in the longitudinal velocity increases as the temperature difference increases. It still exhibits a small-amplitude pulsation oscillation under small temperature differences, while it exhibits a large-amplitude pulsation oscillation under large temperature differences. Similar to transverse velocity, the spectrum of longitudinal velocity has only one fundamental frequency at small temperature differences, and the small amplitude oscillation at this time is a single-periodic oscillation. However, multi-periodic oscillations occur under large temperature differences.

Monitoring the velocity at the high and middle points, it can be observed that the oscillation characteristics of thermo-solutocapillary convection are strongly affected by the temperature difference. As the temperature difference increases (ΔT = 1 − 6), the oscillation amplitude at the monitoring point gradually increases, and the oscillation transforms into large-amplitude pulsation. This indicates that the flow becomes more unstable at higher temperature differences, consistent with previous studies on thermocapillary convection. The analysis of power spectral density (PSD) also reveals interesting information about the oscillation characteristics of capillary convection. Under small temperature differences (ΔT = 1 − 3), the PSD diagram has only one or two main frequencies, indicating the existence of single- or double-periodic oscillations. However, under large temperature differences (ΔT = 6), multiple main frequencies can be observed, resulting in multi-periodic oscillation. This indicates that the flow becomes more complex at higher temperature differences, and multiple frequencies indicate more complex oscillations. In summary, these results indicate that the oscillation characteristics of thermo-solutocapillary convection are strongly affected by temperature difference and become more complex with the increase in temperature difference.

3.2. Effect of Heat Loss on Flow Structure and Interfacial Oscillation Characteristics of Thermo-Solutocapillary Convection

Figure 11a exhibits the contour plot of the internal temperature of the liquid bridge under different heat loss conditions. The results show that the heat loss has little influence on the temperature change in the center vicinity of the system (x = 1.5 to 2.5) via comparing the temperature contours under three heat loss conditions. The temperature profile in this region is almost straight, indicating that the heat transfer pattern near the center is mainly in the form of heat conduction. However, there are some subtle changes in the temperature field of the area due to heat loss in close proximity to the free interface. The temperature near the free interface no longer remains uniform compared to that in the central portion. The non-uniformity increases as the Bi number increases. The closer to the interface, the more apparent the temperature variations. The heat transfer method near the free interface is mainly conducted by convection under the influence of the surface flow, which results in a certain unevenness in the temperature distribution of the region compared with the central part. The heat loss of the liquid bridge magnifies when the Bi number increases. The temperature value at a fixed distance from the interface exhibits a continuously decreasing trend, resulting in a gradual increase in the non-uniformity of temperature distribution near the interface. Figure 11b exhibits the contour plot of the internal concentration under different heat loss conditions. Generally, the concentration value on the interface is smaller than the concentration value at the rest of the same height, and a large concentration gradient is formed at the upper end. Compared with the temperature contour, the concentration contour shows greater heterogeneity. Under the three heat loss conditions, most of the concentration values in the middle part remain almost unchanged. However, the concentration profile near the free interface gradually approaches the upper end of the interface as the Bi number increases, and a large concentration gradient is formed in this region.

Figure 12a exhibits the transverse velocity oscillation distribution of the upper monitoring point “a” under different heat loss conditions, and Figure 12b shows the PSD of transverse velocity. As listed in Figure 12a, the transverse velocity of the upper monitoring point under different heat loss conditions undergoes a stable development process firstly, and then begins to oscillate regularly, and the amplitude gradually increases, eventually entering the high-frequency stable oscillation stage. The yellow area is the selected high-frequency stable oscillation stage, and the arrow refers to an enlarged area map of the area. Figure 12b shows that the monitoring points have only one main frequency under different heat losses, and they all have a single-periodic oscillation. When Bi = 0.2, the average velocity is −0.00227, the average amplitude is 1.2 × 10⁻⁴, and the main frequency is f₁ = 0.053. At Bi = 0.4, the velocity, amplitude, and frequency are −0.00198, 1.1 × 10⁻⁴, and 0.053, respectively. At Bi = 0.6, the velocity, amplitude, and frequency are −0.00174, 8 × 10⁻⁵, and 0.047, respectively. The transverse velocity value and amplitude gradually decrease, while the average oscillation period gradually increases as the Bi number increases. Because the upper monitoring point is located in a region of high temperature and concentration, the fluid flow experiences significant influence from the thermo-solutocapillary convection. Meanwhile, the temperature and concentration of the monitoring point “a” decrease as the degree of heat loss increases (see Figure 11). Therefore, the effects of thermal and solute capillary convection are weakened, resulting in a gradual decrease in flow velocity and amplitude as well as a weakening of oscillation.

Figure 13a exhibits the longitudinal velocity oscillation of the upper monitoring point “a” under different heat loss conditions, and Figure 13b exhibits the PSD of longitudinal velocity. Figure 13a illustrates that the development processes of longitudinal and transverse velocities are different. After experiencing a stable development process, the longitudinal velocity will enter a pulsating oscillation, and the oscillation rule is not as obvious as the transverse velocity. Figure 13b shows that the monitoring point has two main frequencies under different heat losses, both of which are double-periodic oscillations. Additionally, the mean velocity, amplitude, and frequency all decrease as the Bi number increases. Meanwhile, the power of the longitudinal velocity is smaller than the power of the transverse velocity generally, indicating that the transverse velocity changes sharply at the active region of the hot corner.

Figure 14a exhibits the transverse velocity oscillation distribution of the monitoring point “b” under different heat loss conditions, and Figure 14b shows the PSD of transverse velocity. Figure 14a exhibits that the transverse velocity of the monitoring point “b” also first undergoes a stable development process under different heat loss conditions and then begins to oscillate. Meanwhile, the amplitude will gradually increase and eventually enter the high-frequency stable oscillation stage. Figure 14b exhibits that the monitoring points have only one main frequency and are all subject to single-periodic oscillation. Although the oscillation exhibits some pulsating characteristics, its pulsatility is not very strong, and this oscillation can still be considered as a stable oscillation. The average transverse velocity values under different heat loss conditions are 1.24 × 10⁻⁴, 6.85 × 10⁻⁵, and 8.43 × 10⁻⁶, respectively. The amplitudes are 8.8 × 10⁻⁵, 7.1 × 10⁻⁵, and 4.5 × 10⁻⁵, respectively. The oscillation frequencies are 0.053, 0.053, and 0.047, respectively. The transverse velocity value, amplitude, and frequency decrease as the Bi number increases.

Figure 15a shows the longitudinal velocity oscillation of the monitoring point “b” under different heat loss conditions, and Figure 15b shows the PSD of longitudinal velocity. The longitudinal velocity also first undergoes a stable development process, and then begins to oscillate, and the amplitude gradually increases and finally enters the high-frequency stable oscillation phase. The longitudinal velocity value, amplitude, and frequency all decrease as the Bi number increases. Meanwhile, the average oscillation period increases with the increase in the Bi number.

4. Conclusions

This paper presents a numerical investigation into the characteristics of oscillating thermo-solutocapillary convection in toluene/n-hexane mixed solutions. An improved level set method is employed to accurately track every displacement of the interface, while the surface tension is addressed using the CSF model. The study employs the mass conservation level set method to examine the flow and interface oscillation properties under various temperature differences and heat losses, as well as to analyze the temperature, concentration, velocity, and deformation at the interface. The main findings can be summarized as follows:

• In the study of the impact of temperature difference (ΔT = 1 − 6) on oscillatory flows in thermo-solutocapillary convection, the velocity of fluid at the liquid bridge interface increases as the temperature difference rises. Upon reaching a certain threshold (ΔT = 6), the fluid flow at the interface becomes disordered, and the oscillation becomes more intense. The interface deformation under different temperature differences presents an S-type distribution. The velocities at the upper and middle positions of the liquid bridge exhibit different oscillation forms. The transverse velocity at the upper end exhibits a sinusoidal oscillation at a small temperature difference (ΔT = 1 − 3), and this sinusoidal oscillation is a single-periodic oscillation. However, the longitudinal velocity exhibits the pulsating oscillation characteristic with a small amplitude and presents a double-periodic oscillation form. The transverse and longitudinal velocities at the intermediate position exhibit small-amplitude pulsation oscillation, both of which are single-periodic oscillations. The velocity oscillations at the upper and middle positions of the free interface are both pulsating oscillations with large amplitude under the condition of a large temperature disparity (ΔT = 6), and these pulsations are multi-periodic in nature.
• In the investigation on the impact of heat loss on the oscillatory behavior of thermo-solutocapillary convection, the temperature contour distribution within the intermediate region of the liquid bridge (x = 1.5 to 2.5) is relatively uniform, while near the interface, the temperature distribution is no longer uniform. The degree of unevenness increases as the Bi number ranges from 0.2 to 0.6. Furthermore, the concentration values within the concentration field are almost identical in the middle of the liquid bridge. However, the concentration contour near the free interface gradually approaches the upper end as the Bi number increases (Bi = 0.2 − 0.6), thus forming a large concentration gradient at the upper end. In comparison with the adiabatic conditions, the mean velocity, amplitude, and frequency at both the upper and middle positions of the free surface exhibit a gradual reduction as the intensity of heat loss increases. This decrease in dynamic response leads to a weakened oscillation intensity.
• The simulation framework employed in this study makes certain assumptions and simplifications that may not accurately represent real-world conditions. The future research should aim to further explore relevant issues in this field by refining the models used and investigating factors such as different container sizes or shear airflows.