


[Contents]  [Index] 
Cooperation with: K. Böttcher, W. Miller, U. Rehse (Institut für Kristallzüchtung (IKZ) Berlin)
Description: As a result of recent semiconductor technology, there is currently much interest in the growth of highquality indium phosphide (InP) crystals. This is especially true in the fields of optoelectronics and radio frequency electronics. There it is known, for example, that InPbased electronic components can perform much more efficiently than the hitherto favored gallium arsenide (GaAs) or silicon, mainly because InP can operate at higher frequencies.
Since May 2003 we have been engaged in a collaborative project with the Institute for Crystal Growth (IKZ) in Berlin, which has involved comparison of 2D and 3D numerical transientflow simulations performed with different solvers at each institute, as a precursor to experimental growth using the vaporpressurecontrolled Czochralski (VCz) equipment at IKZ (http://www.ikzberlin.de/groups/ag.czhl/pview.php?lang=EN&opt=0&no=07). This method of growth is a variant on the liquidencapsulated Czochralski (LEC) growth method; an important additional feature of the VCz method is that the crystal, melt, and encapsulant are placed in an insulated inner chamber, which can significantly reduce potentially damaging thermal stresses within the grown crystal. We have previously applied the method to modeling the growth of GaAs crystals, as documented in last year's report [1]. Subsequent results can be found in [2], [3]. A similar model to the one described therein can be applied to InP growth, i.e. the liquid mechanics are again described by a coupled thermohydrodynamic system comprising the incompressible NavierStokes equations and the heattransport equation. A fundamental difference applies to the boundary conditions for temperature, specifically on the encapsulant/melt interface and the crucible wall . For the InP study, we have applied a Dirichlet condition on both boundary sections, based on advectionfree global simulation data supplied by IKZ. This contrasts with the more simplified conditions of adiabaticity on and Dirichlet value dependence derived from simple functions on , as assumed for the GaAs study.

In our mathematical model, the flow variables are nondimensionalized in a standard fashion, mainly using the crucible radius R_{C} (as characteristic length), the melt kinematic viscosity , the melt density , the temperature difference T between the crystal/melt interface and the rim of the crucible, and the angular velocities of the crystal and crucible ( and , in turn). Industryrelevant values of these quantities are shown in Table 1 above, as well as the Prandtl number Pr and thermal Rossby (or Richardson) number Ro. Finally, a nondimensional crystal radius r_{X} = 0.5 and a nondimensional height h_{T} = 0.6 were chosen. As with the GaAs project, we have performed computations using both cylindrical and realistic geometry. Also our initial efforts have been directed towards the case of isorotation, with = , and this holds for the examples shown below. The Reynolds number is appropriately defined as: Re = R_{C}^{2}/.
First, in the case of cylindrical geometry, we have performed axisymmetric, and more recently, threedimensional (3D) meltflow simulations for Reynolds numbers ranging from 1000 to 8000, all with Ro = 1 , Pr = 0.015. Although these values are typically well below those required by industry (whose values lie inside the turbulent regime), the rich structure in the resultant flow behavior seems worthy of analysis in its own right. Moreover, it provides a good basis for comparing results obtained with our finite element solver NAVIER, and those provided by STHAMAS3D, a finitevolume code with upwinding implemented at IKZ. The latter is more dissipative than ours and may be better suited to computing flow properties more relevant to industry, especially timeaveraged quantities such as crystal/melt interface shape, and mean frequencies. On the other hand, NAVIER seems to be better suited to capturing transient effects accurately. As in the case for the GaAs research, the steady axisymmetric state was found to be most unstable to 3D oscillatory instabilities; for isorotation with = for example, we have found that, while the 2D solver first exhibits transient behavior for Re 4500 (corresponding to 1.18 rpm), in the 3D simulation this value is just below 3000 (around 0.79 rpm).
In fig1_vcz.eps, we can see the time history plots for a selection of point data, as well as the Nusselt number on the crystal/melt interface and the mean difference between successive solutions; the relatively long ``swingin'' phase which precedes oscillatory motion may suggest nearness to the critical Reynolds number. It is also noticeable from these plots that the high frequencies at this stage are purely associated with the threedimensionality, with twodimensional frequencies being very short, in comparison.
fig2_vcz.eps depicts the timeaveraged flow and temperature profiles in an (almost) arbitrary vertical plane passing through the crucible and containing the axis. Here the time averaging was based on using a nondimensionalized end time of 200 (or equivalently, running for 31.8 revolutions in dimensional terms). In fact, this solution is very similar to the axisymmetric component of the flow and reinforces the point regarding the relatively weak threedimensional presence. Examination of the Fourier angle modes for various flow quantities reveals that the mode with azimuthal wavenumber 4 is the dominant 3D mode, although this is virtually indiscernible from crosssectional slices of the crucible, so dominant is the axisymmetry. Another conclusion is that the flow is dominantly rotational with buoyancy playing a minor rôle here; the strongest nonrotational influence appears to stem from the large temperature gradient near the crystal/melt/encapsulant triple line (annulus).
We have also recently applied our axisymmetric solver to simulating melt flows in crucibles of a more realistic design. The geometry here principally consists of a cylindrical flattopped rim section joined to a spherical section, which itself is joined to a spherical cap of relatively large radius. Although similar qualitative features are found in comparison with using cylindrical geometry, early indications suggest that the former type is more destabilizing. In fig3_vcz.eps, for example, the solution is presented for Pr = 0.015 , Ro = 1, and Re = 3500 (0.92 rpm, approx.) portraying the timeaveraged temperature and projectedvelocity fields, based on a dimensionless end time of 100 (15.9 revolutions).
In conclusion then, we have computed axisymmetric and 3D steady/transient
meltflow solutions for crucibles with cylindrical geometry, and axisymmetric
steady/transient solutions using realistic geometry.
In general, all of the results tend to suggest that, for Reynolds numbers up to
8000 at least, the flow is dominantly rotational, with the strongest
nonrotational influence stemming from the triple line.
For future work, we are principally interested in
simulating transient 3D solutions using realistic geometry
(in parallel with IKZ). Before this however, a realistic design of the meniscus
at the triple line will be added, and suitable 3D meshes constructed.
References:



[Contents]  [Index] 