Journal of Crystal Growth 385 (2014) 61–65
Contents lists available at ScienceDirect
Journal of Crystal Growth journal homepage: www.elsevier.com/locate/jcrysgro
The relative contributions of thermo-solutal Marangoni convections on flow patterns in a liquid bridge H. Minakuchi a,n, Y. Takagi b, Y. Okano b, S. Gima a, S. Dost c a
Department of Mechanical Systems Engineering, University of the Ryukyus, 1 Senbaru, Nishihara, Okinawa 903-0213, Japan Department of Materials Engineering Science, Graduate School of Engineering Science, Osaka University, 1-3 Machikaneyama, Toyonaka, Osaka 560-8531, Japan c Crystal Growth Laboratory, University of Victoria, Victoria, BC, Canada V8W 3P6 b
art ic l e i nf o
a b s t r a c t
Available online 14 June 2013
A numerical simulation study was carried out to investigate the relative contributions of thermal and solutal Marangoni convections on transport structures in a liquid bridge under zero gravity. The liquid bridge in the model represents a three dimensional half-zone configuration of the Floating Zone (FZ) growth system. Three dimensional field equations of the liquid zone, i.e. continuity, momentum, energy, and diffusion equations, were solved by the PISO algorithm. Computations were performed using the open source software OpenFOAM. The numerical simulation results show that the flow field becomes three-dimensional and time-depended when the solutal Marangoni number is larger than the critical value. It was also shown that not only flow patterns but also the azimuthal wave number (m) changes due to the competing contributions of thermal and solutal Marangoni convective flows. Crown Copyright & 2013 Published by Elsevier B.V. All rights reserved.
Keywords: A1. Computer simulation A2. Floating zone technique Half-zone Marangoni convection
1. Introduction When there is a free surface in a melt, the variation in the surface tension drives a convective flow along the free surface in the melt. This flow is known as Marangoni convection, and may be produced by temperature (thermal Marangoni convection) and/or concentration gradients (solutal Marangoni convection) in the melt. It is known that such Marangoni convective flows may become unstable much easier than the buoyancy driven natural convection. When a free surface exists in a crystal growth technique such as the Floating-zone (FZ) technique, unstable Marangoni flows that develop in the melt (liquid bridge) may induce undesirable growth striations in the grown crystals. FZ crystal growth being a containerless technique (in which the melt is solidified through a liquid bridge with a free surface) has the advantage of producing crystals free from crucible contamination. However, when larger liquid bridges are selected in order to grow large crystals, or/and when temperature and concentration gradients are steeper, Marangoni convection may become stronger and lead to flow instabilities in the melt. In the growth of SixGe1−x single crystals, a concentration gradient in the melt leads to further segregation. Such a concentration gradient gives rise to a strong solutal Marangoni convection in the melt compared with thermal Marangoni convection. It is important to know the relative contributions of these convections in the
n
Corresponding author. Tel.: +81 98 895 8619; fax: +81 98 895 8636. E-mail address:
[email protected] (H. Minakuchi).
growth of SixGe1−x crystals by FZ [1,2]. From this point of view, the present study investigates the role of thermal and solutal Marangoni convections in the Floating Zone system numerically. Most modeling studies in literature have considered thermal Marangoni convection only [3–7]. Those considering solutal Marangoni convection were few [8], including an experimental work on the Si–Ge crystal growth that has reported the adverse effect of solutal Marangoni convection in crystal growth [9]. To this end, a numerical simulation study would be very handy to shed light on the relative contributions of such convective flows. In this study, this is done through three-dimensional numerical simulations for the thermal and solutal Marangoni convections in the melt of the Si–Ge FZ growth system, and their effects on transport structures of the melt are examined. In order to make the relative contributions of thermo-solutal Marangoni convections distinguishable, and not to involve the contribution of natural convection, the gravity level was taken zero (g¼0) in the model.
2. Mathematical model In FZ growth of single crystals, various factors such as growth rate, interface deformation, the concentration condition at the liquid–solid interface and the temperature condition at the free interface may play roles as described in [10]. However, since the objective of the present study is to shed light on the relative contributions of thermo-solutal Marangoni convection on flow structures of the melt under zero gravity, we have simplified the
0022-0248/$ - see front matter Crown Copyright & 2013 Published by Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.jcrysgro.2013.05.042
62
H. Minakuchi et al. / Journal of Crystal Growth 385 (2014) 61–65
ΔC ΔT θ μ ν ρ s
Nomenclature a Asp C D L m Ma N p Pr r Sc t T ! V z α θ
radius of liquid bridge [m] aspect ratio (¼L/a) molar fraction of Si [–] diffusion coefficient of SixGe1−x [m2/s] length of free surface [m] azimuthal wave number Marangoni number mesh number pressure [Pa] Prandtl number radial coordinate [m] Schmidt number time [s] temperature [K] velocity vector [m/s] axial coordinate [m] thermal diffusivity [m2/s] azimuthal coordinate [rad]
Subscripts 0 c C h p T R Z θ
model as much as possible by not taking into account some of these factors to remove their possible influences. With this in mind, we have made the following assumptions in the model: (i) the Si–Ge melt is an incompressible, Newtonian fluid mixture, (ii) the solid/liquid interfaces are flat, and (iii) the system is under zero gravity. Since the gravity is taken zero, the liquid bridge remains cylindrical (and thus the liquid/gas interface is flat), and natural convection does not develop in the melt (due to the absence of the gravitational body force). The three-dimensional half-zone of the FZ liquid bridge is considered as the liquid phase of the model domain (its schematic description is given in Fig. 1). Under these assumptions, the following governing equations of the liquid phase are obtained respectively from the overall mass conservation, the balance of momentum, the balance of energy and the conservation of mass of species: Continuity ! ∇⋅ V ¼ 0
ð1Þ
Momentum ! ! ! ! ∂V 1 þ V ⋅∇ V ¼ − ∇p þ νΔ V ∂t ρ
initial cold solutal hot period thermal R-direction Z-direction θ-direction
Upper-disc (z¼ L): Along the upper boundary we consider the non-slip condition on the flow velocity field, and assumed that temperature is prescribed at T¼Tc. The melt concentration is assumed as the pure silicon at this boundary: C¼1
ð5Þ
Lower-disc (z¼0): Along the lower boundary, similarly we also use the non-slip condition on the flow velocity field, and assumed that temperature is prescribed at T ¼Th. The melt concentration on this boundary is taken as the pure germanium C¼0
ð6Þ
The concentration boundary conditions imply that we assume the melt is in equilibrium with the discs (solids) concentrations during melting and growth. Free surface (r ¼a): Along the free surface we use the following conditions. The contribution of Marangoni convection was taken into account using the following tangential force balance at the free surface: Vr ¼ 0
ð2Þ
Energy ∂T ! þ V ⋅∇T ¼ αΔT ∂t
concentration difference [–] temperature difference [K] azimuthal coordinate [rad] viscosity [Pa s] kinematic viscosity [m2/s] density [kg/m3] surface tension [N/m]
ð7Þ
μ
∂V z ∂s ∂T ∂s ∂C þ ¼− ∂T ∂z ∂C ∂z ∂r
ð8Þ
ð3Þ
∂ Vθ 1 ∂s ∂T ∂s ∂C μ r þ ¼− ∂r r r ∂T ∂θ ∂C ∂θ
ð4Þ
Table 1 Parameters and properties of the SixGe1−x system [13].
ð9Þ
Mass transport ∂C ! þ V ⋅∇C ¼ DΔC ∂t
Boundary conditions: We consider the following boundary conditions in the computational domain.
Fig. 1. Schematics of the liquid bridge.
Property
Symbol
Value
Temperature coefficient of surface tension Concentration coefficient of surface tension Radius of liquid bridge Length of liquid bridge Aspect ratio Thermal diffusivity Kinematic viscosity Diffusion coefficient Prandtl number Schmidt number
∂s/∂T [N/m K]
(0,−2.20,−5.51) 10−6
∂s/∂C [N/m]
0−2.57 10−5
a [m] L [m] Asp α [m2/s] ν [m2/s] D [m2/s] Pr Sc
1.00 10−2 5.00 10−3 0.5 2.20 10−5 1.40 10−7 1.00 10−8 6.37 10−3 14.0
H. Minakuchi et al. / Journal of Crystal Growth 385 (2014) 61–65
The radial temperature and concentration gradients are zero along the free surface: ∂T ∂C ¼ ¼0 ∂r ∂r
ð10Þ
The associated dimensionless thermal and solutal Marangoni numbers are defined as ∂s ΔTL ð11Þ MaT ¼ − ∂T μν MaC ¼
∂s ΔCL ∂C μν
ð12Þ
63
The initial conditions were: no flow (zero flow velocity) and prescribed constant temperature and concentration values (T¼ Tc and C¼1). Eqs. (1–4) were discretized by the Finite Volume Method on a collocated grid system, and solved by the PISO algorithm. The computation was carried out using the open source software OpenFOAM. Euler scheme, QUICK scheme and Gauss linear scheme were applied respectively to the terms involving time derivative, divergence, and Laplacian in the governing equations. Physical properties of the SixGe1−x system and their values are listed in Table 1. Computations were performed using the NR Ny NZ ¼ 30 120 40 mesh. Non-uniform grids in the region near the upper and lower discs and free surface were used as shown in Fig. 2. Details, such as the model description, numerical scheme, and validity of simulation, can be found in [11,12].
3. Results and discussion 3.1. Thermal Marangoni convection
Fig. 2. Computational grid.
Figs. 3 and 4 show the results when the thermal Marangoni convection is considered alone (at MaT ¼1072 or 2679 and MaC ¼0). In these figures, the temperature distribution, the SixGe1−x concentration distribution, and the flow fields in the horizontal (z¼L/2) and vertical
Fig. 3. Computed temperature, melt concentration, and flow intensity values at MaT ¼1072 and MaC ¼ 0.
Fig. 4. Computed temperature, melt concentration, and flow intensity values at MaT ¼ 2679 and MaC ¼ 0.
64
H. Minakuchi et al. / Journal of Crystal Growth 385 (2014) 61–65
planes are presented. Shear instability gave rise to the transition of an axisymmetric flow to a 3-D steady flow with the azimuthal wave number of m¼3, as shown in these figures. However, as can be seen from Fig. 4, the temperature field is still axisymmetric (almost), while the concentration and flow fields exhibit a three-dimensional behavior. This is because the Schmidt number is much larger than the Prandtl number.
3.2. Solutal Marangoni convection Figs. 5 and 6 present the computed melt concentration distribution in the horizontal and vertical planes when only the solutal Marangoni convection is considered (at MaT ¼0 and MaC ¼ 1072 or 1107). At a smaller solutal Marangoni number than the critical number (about MaC ¼ 1107), the flow pattern is still axisymmetric. On the other hand, at a larger solutal Marangoni number than the critical number, the flow becomes unstable, and a 3-D pulsating oscillation of m¼4 appears as shown in Fig. 6. The concentration field oscillates at a period of 64 s. 3.3. Thermo-solutal Marangoni convection
Fig. 5. Melt concentration distribution at MaT ¼ 0 and MaC ¼1072.
Fig. 7 presents the computed results for the case of thermosolutal Marangoni convection (taking into account both thermal and solutal Marangoni convections at MaT ¼2679 and MaC ¼1107 or 1250). As shown in Fig. 7(a), although the solutal Marangoni number (MaC ¼ 1107) was larger than the critical number of the case of only-solutal Marangoni convection, the flow pattern is almost axisymmetric. On the other hand, at MaT ¼2679 and MaC ¼1250, unlike the flow patterns of the individual cases of only-thermal and only-solutal Marangoni convections, the flow pattern in this
Fig. 6. Time dependency of the melt concentration at MaT ¼0 and MaC ¼ 1107 (tp ¼64 s).
Fig. 7. Time dependency of the melt concentration at MaT ¼ 2679 and MaC ¼ 1107 (a) and 1250 (b) (tp ¼ 16 s).
H. Minakuchi et al. / Journal of Crystal Growth 385 (2014) 61–65
65
Fig. 8. Time dependency of the flow velocity (left) and the melt concentration (right) at MaT ¼ 0 or 2679 and MaC ¼ 0 or 1107.
Fig. 9. Time dependency of the melt concentration at MaT ¼ 2679 and MaC ¼ 982, 1107 and 1250.
combined case becomes rotating oscillatory with m¼7 as shown in Fig. 7(b). Due to the competing contributions of the thermal Marangoni and solutal Marangoni convections, not only the flow patterns but also the azimuthal wave number (m) change. Fig. 8 shows the time dependency of the flow strength (left) and the melt concentration (right) at the point of (R, y, Z)¼(0.999a, 0, 0.5L). As shown in Fig. 8 (left), the flow velocity magnitude in the combined case (thermo-solutal Marangoni convection) is increased compared with those of the individual cases, and it appears that it has reached a steady value. However, because of low amplitude and high frequency vibrations as shown in Fig. 8 (right), the flow pattern remained almost axisymmetric as can be seen from Fig. 7(a). Fig. 9 shows the time dependency of the melt concentration at MaT ¼2679 and MaC ¼982, 1107 and 1250 at the point of (R, y, Z)¼ (0.999a, 0, 0.5L). In the case of smaller solutal Marangoni numbers than MaC ¼982, the convective flow was steady and axisymmetric. At higher solutal Marangoni numbers (MaC ¼1107 and 1214), the convection became unstable. However, the flow pattern was still almost axisymmetric as can be seen from Fig. 7(a). At a larger solutal Marangoni number than the critical number (MaC ¼ 1250), the convection became a 3-D oscillatory at a period of 16 seconds. 4. Conclusions The numerical simulation study carried out to examine the role of Marangoni convection in a three-dimensional half-zone liquid bridge led to the following conclusions: 1. The convection at MaT ¼0 and MaC ¼1107 becomes unstable, and a 3-D pulsating oscillation of m ¼4 occurs. 2. Unlike the flow patterns in the individual cases of thermal and solutal Marangoni convections, the flow pattern of the thermo-
solutal Marangoni convection case becomes oscillatory (rotating oscillatory of m ¼7). 3. The critical solutal Marangoni number of thermo-solutal Marangoni convection is larger than that of only-solutal Marangoni convection.
Acknowledgments This work was partially supported by a Grant-in-Aid for Scientific Research (B) (no23360343) from the Ministry of Education, Culture, Sports, Science and Technology of Japan. References [1] H. Minakuchi, Y. Okano, S. Dost, Journal of Crystal Growth 266 (2004) 140. [2] T.P. Lyubimova, R.V. Skuridin, I.S. Faizrakhmanova, European Physical Journal Special Topics 192 (2011) 41. [3] C.W. Lan, J.H. Chian, Journal of Crystal Growth 230 (2001) 172. [4] Z. Zeng, H. Mizuseki, K. Simamura, T. Fukuda, K. Higashino, Y. Kawazoe, International Journal of Heat and Mass Transfer 44 (2001) 3765. [5] M. Lappa, R. Savino, R. Monti, International Journal of Heat and Fluid Flow 44 (2001) 1983. [6] H. Minakuchi, Y. Okano, S. Dost, International Journal of Materials and Product Technology 22 (2005) 151. [7] K. Li, B. Xun, N. Imaishi, S. Yoda, W.R. Hu, International Journal of Heat and Fluid Flow 29 (2008) 1190. [8] L.M. Witkowski, J.S. Walker, Physics of Fluids 14 (2002) 2647. [9] T.A. Campbell, M. Schweizer, P. Dold, A. Crӧll, K.W. Benz, Journal of Crystal Growth 226 (2001) 231. [10] T.P. Lyubimova, R.V. Skuridin, I.S. Faizrakhmanova, Journal of Crystal Growth 303 (2007) 274. [11] H. Minakuchi, Y. Takagi, Y. Okano, K. Mizoguchi, S. Gima, S. Dost, Journal of Advanced Research in Physics 3 (2012) 011201. [12] H. Minakuchi, Y. Takagi, Y. Okano, T. Nosoko, S. Gima, S. Dost, Transactions of JSASS, 10, Aerospace Technology, Japan, 2012, Ph_15. [13] S. Abbasoglu, I. Sezai, International Journal of Thermal Sciences 46 (2007) 561.