Axial Vibration Control Technique for Crystal Growth from the Melt: Analysis of Vibrational Flows’ Behavior

Abstract

A problem of efficacy of crystal growth methods for crystallization from solutions or melt has been investigated. The axial vibrational control (AVC) technique was considered as a perspective method to manage both heat-mass transfer and chemical component composition of the melts in the case of crystallization of complex chemical compounds. Numerical modeling and the search for generalized dependencies made it possible to predict the AVC parameters that provide optimal heat and mass transfer modes for creating flat liquid-solid interfaces, as well as the component composition of dissociated melts of various chemical compounds—Ge, NaNO₃, CdTe.

1. Introduction

Over the past 20 years, in the field of functional materials, there has been a market shift towards thin-film and nano-sized structures [1]. Nevertheless, the market for single crystals still continues to grow [2]. In spite of the progress in polycrystal technology, many properties which are found in single crystals cannot be replicated in polycrystals. In many practical applications, only single crystal materials are able to provide the required and record-breaking functional characteristics [3,4].

To date, single crystals of semiconducting materials are one of the most widely researched and used materials. In addition to Si and Ge, the single crystals of III–V semiconductors, due to their outstanding optical and electronic properties, are an integral part of devices for application in fiber-optic communication, wireless and satellite communication, solid-state lighting, etc. [5,6].

In this regard, the tasks of both improving quality, primarily structural, and reducing the cost of crystals are fundamental. Taking into account the fact that more than 90% of all mass-production single crystals are grown from the liquid phase (solution or melt growth methods), the issue of improving these methods is extremely relevant.

Numerous studies on improving of crystal growth methods by various types of mixing (ultrasonic [7,8], ACRT—accelerated crucible rotation technique [9], CVS—coupled vibrational stirring [10], rotating baffle [11], HFR—heat field rotation [12], TMF—travelling magnetic field [13], AHP—axial heat flux technique [14]) were mainly aimed at effective control of heat and mass transfer, creating conditions for obtaining a flat interface in the case of directional crystallization methods.

At the same time, the issue of control of chemical thermodynamics states in growth systems has been considered extremely infrequently [15].

A study of the influence of low-frequency vibrations on the equilibrium of chemical reactions [16] showed that a clear effect of chemical equilibrium control was possible. In ref. [17], when analyzing the crystals growth processes of chemical compounds, an explanation for this effect was given. It has been established that when an AVC is introduced into a liquid through harmonic vibrations of a solid chemically inert flat cylindrical disk, which has a sharp edge on a cylindrical generatrix, such high-energy flows are generated in the region of the sharp edge that the energy of viscous dissipation of liquid layers becomes comparable to the chemical interaction energy of clusters in the melt. The transformation of the mechanical energy of a moving liquid into chemical energy leads to a change in the chemical nature of the liquid phase, in particular, to a change in the component composition and destruction of associates in the melt. From this point of view, the AVC technique can be considered not only as an effective tool for organizing the required heat and mass transfer, but also as a tool for controlling the thermodynamic properties of the liquid phase.

The efficacy of the AVC technique application for growth of both dielectric [18,19], and semiconductor crystals [20] has been demonstrated. The application of the AVC technique made it possible to both improve the structural perfection of the grown crystals and increase the volumetric growth rate in several times. However, for each specific crystal and each method, the authors were forced to experimentally optimize the AVC modes to achieve the desired result.

Taking into account the above, it would be extremely desirable to have a generalized dependence of the influence of the parameters of the AVC method on the properties of flows in the liquid phase and their component composition. Optimization of the AVC method for growing crystals is a multi-parameter problem in which it is possible to distinguish the control parameters of the AVC technique (oscillation amplitude and frequency, disk geometry) and the parameters of a specific technology (properties of the crystallizing melt, geometry of the crystal growing method). To achieve the generalized dependence, first of all, it is necessary to establish the dependence of the speeds of vibration flows on the amplitude, frequency, density, and viscosity of the liquid for a cylindrical vibrating body with a different radius of edge rounding. To achieve this, this study excluded all factors such as free sur-face, thermal convection, and crystallization front, which are always present in a real environment of growth process and interfere with an objective assessment of the characteristics of vibration flows.

2. Numerical Model

To understand the nature of the vibratory flow formation and properties, we have excluded from consideration the thermal convection which always accompanies the melt crystal growth process by conducting the system analysis in ideal isothermal conditions. Numerical simulation was conducted using the ANSYS Fluent package and germanium. Cadmium telluride and sodium nitrate melt were considered as model substances. The AVC technique was realized by low-frequency harmonic oscillations of a cylindrical disk submerged into the melt. The disk is specified by three critical parameters: (1) diameter—D, (2) thickness—d, and (3) radius of the edge curvature between a cylindrical generatrix and a faceplate—R. It had been demonstrated that the R-value sufficiently influenced on average flows’ generation [18].

2.1. Model of Fluid Motion under the AVC Action

The analysis of the AVC technique was carried out by numerical modeling in a two-dimensional axisymmetric approximation using the geometric model shown in Figure 1. The grid model (see Figure 1 right) contains 187,000 cells for modeling half of the axially symmetric model with an average size of 0.3 mm in the region of the liquid phase with compaction of cell layers at the boundaries. The first mesh layer at the disk surface is 1.7 μm thick and the minimum cell size on the disk rounding is 4.9 μm with a gradual expansion into the main design volume. The dimensions of the geometric model are presented in

Table 1. Geometrical parameters of the numerical model.

Parameter	Value, mm
Simulating area length	100
Simulating area height	100
Disk diameter, D	80
Disk thickness, d	10
Edge radius, R	0.1; 0.2; 0.4; 0.8
Distance from the disk to the bottom and top walls	45
Disk oscillation amplitude, A	0.080–0.500
Frequency oscillation, f	10–30 Hz

.

The simulating area was divided into five zones (Figure 1 center). In a zone surrounding the vibrator—disk zone—the grid was moving with the vibrator as well. In the two zones belonging to the upper and lower ends of the crucible—stationary zones—grids were non-moving. There were two zones between the moving disk zone and the non-moving stationary zones at the crucible ends. In these zones, the grid was stretched and compressed with disk motion.

The velocity and pressure fields were calculated using the Navier–Stokes equation, supplemented by the continuity equation, at constant density and viscosity:

$$\rho \left(\frac{\partial \overrightarrow{v}}{\partial t}+\left(\overrightarrow{v}·\nabla \right)\overrightarrow{v}\right)=-\nabla p+\mu ∆\overrightarrow{v};$$

(1)

$$\nabla ·\overrightarrow{v}=0,$$

(2)

here $\rho$ represents the density of the melt, $\overrightarrow{v}$ represents the velocity of the moving melt, $p$ represents the melt pressure, $\mu$ represents the melt dynamic viscosity coefficient.

The model was calculated under isothermal conditions:

$$T=\mathrm{const},$$

(3)

here $T$ represents the melt temperature.

Boundary conditions for a stationary wall were specified using no-slip boundary condition:

$${\overrightarrow{v}}_n=0,$$

(4)

here $\overrightarrow{n}$ represents the normal direction of the surface.

The oscillation of the disk along the z-axis was specified using the boundary condition for changing the speed:

$$\left\{\begin{array}{c}{\mathit{v}}_z=\mathit{A}\mathit{\omega }·\mathit{cos}\left(\mathit{\omega }t\right);\\ {\mathit{v}}_r=0,\end{array}\right.$$

(5)

$$\mathit{\omega }=2\pi \mathit{f},$$

(6)

here ${\mathit{v}}_z$ represents the velocity component of oscillating disk along z-axis, ${\mathit{v}}_r$ represents the velocity component of oscillating disk along r-axis, $\mathit{A}$ represents the vibration amplitude, $\mathit{\omega }$ represents the angular frequency and $\mathit{f}$ represents the vibration frequency.

The numerical integration of the system of differential equations was carried out using the finite volume method. The search for a solution was carried out using the iterative Gauss–Seidel method using the pressure-velocity separated SIMPLE scheme. The gradients were discretized using the least squares method. Spatial discretization of velocity and pressure was carried out according to a second-order scheme. Time discretization was carried out according to a first-order scheme. The disk movement was accompanied by the use of a dynamic mesh.

In this work, we used a UDF file that describes the motion of the mesh according to the following iterative formula:

$$x_j^i=x_j^{i-1}+v_z^i·\frac{x_j^{i-1}-x_0^{i-1}}{x_n^{i-1}-x_0^{i-1}}·∆t,$$

(7)

here $x_j^i$ represents the position of node j (from 0 to n) at the current time step i; $x_j^{i-1}$ represents the position of node j (from 0 to n) at the previous time step i − 1; $x_0^{i-1}$ represents the position of the fixed boundary of the deformable mesh zone at the previous time step i − 1; $x_n^{i-1}$ represents the position of the moving boundary of the deformable zone at the previous time step i − 1; $v_z^i$ represents the velocity of the disk from Equation (5) at the current time step I; $∆t$ represents the time step.

The entire deformation procedure is determined by a single value of the current displacement of the disk relative to the average position; old and new coordinates of individual grid nodes are not stored separately anywhere; they are simply recalculated using a formula.

2.2. Material Properties

NaNO₃, Ge and CdTe melts were used as model liquid phases. Their physicochemical properties are presented in

Table 2. Properties of the melts under study at T = T_m + 5 K.

Property	Value
	NaNO3	Ge	CdTe
Melting temperature, Tm (K)	579.95	1211.4	1369
Density, ρ (kg/m3)	1903	5490	5661.8
Dynamic viscosity, μ (Pa s)	2.2 × 10−3	7.14 × 10−4	2.33 × 10−3

.

3. Results

3.1. AVC Flows Velocities

For the selected materials, non-stationary calculations were carried out with a time step equal to 1/400 of the oscillation period until a stationary regime of melt flow was established, that is, until the velocity values averaged over the oscillation period reached a constant value. Distributions of the flow velocity of melts of Ge, CdTe, and NaNO₃ were obtained at various values of the amplitude (A) and frequency (f) of disk oscillation, as well as at various radii of the disk edge (R) (

Table 3. Modules of the maximum melt velocity along the central axis at various values of amplitude, frequency, disk edge radius and intensity.

f, Hz	A, μm	R, mm	I, kg0.5 m1.5 s−1	v, mm·s−1
Ge
10	125	0.1	8.00	0.126
10	200	0.1	32.77	2.240
20	100	0.1	16.39	0.744
20	110	0.1	21.81	1.199
25	80	0.1	13.11	0.535
25	100	0.1	25.61	1.370
25	110	0.1	34.08	2.217
25	125	0.1	50.01	3.635
30	110	0.1	49.08	3.376
25	100	0.2	12.80	0.587
25	100	0.4	6.40	0.096
25	100	0.8	3.20	0.060
25	50	0.05	6.40	0.096
CdTe
10	125	0.1	2.49	0.064
10	200	0.1	10.20	0.657
10	300	0.1	34.42	2.683
20	100	0.1	5.10	0.217
25	80	0.1	4.08	0.091
25	100	0.1	7.97	0.474
25	125	0.1	15.56	1.114
25	500	0.1	996.03	24.86
25	500	0.2	498.01	20.18
25	500	0.4	249.01	14.55
25	500	0.8	124.50	8.757
NaNO3
25	500	0.1	611.57	21.67
25	500	0.2	305.79	16.27
25	500	0.4	152.89	9.655
25	500	0.8	76.45	4.857

).

At the beginning of oscillations (over several periods of oscillation), the disk forms a solenoidal melt flow field in space and sinusoidal in time (Figure 2), which we normalized for clarity to the maximum speed of oscillation of the disk Aω. It can be seen that near the edge of the disk the speed of movement of the melt increases by approximately 10 times the speed of movement of the disk, and in the area between the disk and the side wall the speed of movement of the melt is opposite to the movement of the disk. However, over time, due to the viscosity of the melt, a thin layer of liquid should appear on the side of the disk, which will move in the same direction as the movement of the disk, which will lead to a break in the existing streamlines and the formation of a controlled vortex flow of the melt.

The closest zone to the surface of the disk (Figure 3) forms a thin viscous layer of melt, which moves in the same direction as the movement of the disk. When reaching the edge of the disk, this layer is divided into two components. One part of the melt layer, according to the oscillatory motion of the disk (Figure 2), changes the direction of movement to the opposite and partially closes near the edge of the disk, forming a primary vortex. The second part of the melt layer, acquiring a fairly high speed at the edge of the disk, is directed at an angle towards the flow around the surface of the disk, flowing into the main mass of the melt, causing, together with the oscillatory component of the field, secondary vortices (Figure 4). There are two such vortices—between the disk and the side surface of the chamber (small vortex) and a vortex in the main region of the melt, which, in the case of growing crystals using the technique of axial low-frequency vibrations, reaches the crystallization front. Therefore, we are most interested in correlating the flow rate of the secondary vortex with the vibrational parameters and physical properties of the melt.

Analysis of fluid flows in normalized scale v/(Aω) at the edge of the disk showed that in the case of more viscous CdTe the thickness of the boundary layer on the side surface of the disk is significantly greater comparing to Ge melt (Figure 3). However, the scale of the edge vortex is almost the same in both cases.

It is also clear that for less viscous germanium, inertial effects are much stronger. Hence, the higher speeds of jet flows from the corner, which create circulation in the main volume of the melt.

The mechanism of vortex generation is the tear off of flows when flowing around an edge of large curvature (small radius of curvature). Due to the tear off of flows during the forward and reverse disk motion, the liquid in the edge region moves along different trajectories, resulting in the formation of an edge vortex. The scale of the vortex is determined by the region of separated flow at the disk edge.

To analyze the influence of oscillatory parameters and physical properties of the substance on the nature of the flow, the modulus of the maximum flow velocity along the central axis was taken as a flow characteristic. We tried to find a power-law relationship (n_i) between the melt velocity and amplitude A, angular frequency ω, disk edge radius R, melt density ρ and dynamic viscosity µ:

$$v ~ {\mathit{A}}^{n_1}{\mathit{\omega }}^{n_2}{\mathit{R}}^{n_3}{\rho }^{n_4}{\mu }^{n_5}.$$

(8)

In the case of laminar stationary motion, the data must fit into a linear relationship. As a selecting result, the degrees of n_i using regression analysis, we revealed (Figure 5) that the speed of the melt motion should be proportional to the expression, which for brevity we will denote as I and will call the intensity of oscillations:

$$v ~ \frac{{\rho }^{0.5}{\mathit{A}}^3{\mathit{\omega }}^2}{\mathit{R}\mu } ~ \mathit{I}.$$

(9)

It can be noted that the full set of points we calculated (Figure 5) at I intensity is divided into a linear region at low values of the amplitude and frequency of oscillations and a logarithmic region at high values amplitude and frequency of oscillations. The linear region corresponds to the laminar flow regime without the vortex shedding from the edge of the disk. The logarithmic region coincides with the beginning of the transition to the turbulent regime. At the edge of the disk, one can observe a loss of flow stability, which leads to the development of small vortices. These vortices, with further development, form a vortex street, which transport into a secondary vortex under the disk.

In the region of low frequencies and vibrations (Figure 6), in which the points of Ge and CdTe fall, a linear dependence is observed with a general slope of the straight line. The role of the shift of lines will be discussed in Section 4 when analyzing the obtained dependence. We observe that the linear dependence, hence the laminar flow regime, remains in the range of oscillation intensity from 10 to 50 kg^−0.5 m^1.5 s⁻¹.

At high frequencies and amplitudes of oscillations (Figure 5), which included points of CdTe and NaNO₃ at f = 25 Hz and A = 500 μm when R changed from 0.1 mm to 0.8 mm, it is clear that the melt flow regime begins to become turbulent. This is explained by the appearance of an uneven velocity field in the secondary vortex, as well as an increase in the dissipation energy in the area under the disk (Figure 7 and Figure 8). Separation vortices begin to form here, which are transferred from the edge of the disk to the region of the secondary vortex, leading to a change in the speed of the secondary vortex. This is most clearly visible in the velocity field of CdTe melt at R = 0.1 mm, where several small vortices (blue dots) can be seen that leave the disk edge and propagate into the central flow. Also, in the central region under the disk, a reverse flow occurs more actively, which also leads to a redistribution of velocity on the axis. It can be seen that with a decrease in the radius of the disk edge, this regime stabilizes, gradually reducing to a laminar flow regime and more potential type of the flow.

3.2. AVC Thermodynamic Control

We carried out a detail analysis of R influence on generation and on parameters of steady vibrational flows as well as on viscous dissipation energy (P_w), determined by the expression:

$${\mathit{P}}_w=\mu {{\sigma }_{ij}}^2;$$

(10)

$${\sigma }_{ij}=\nabla \overrightarrow{v}+{\left(\nabla \overrightarrow{v}\right)}^T.$$

(11)

During crystal growth processes, the average vibrational flows are responsible for heat-mass-transfer, while P_w determines structure characteristics of the melt of a complex chemical composition [17]. Numerical simulation showed that the maximum (averaged over the oscillation period) velocity of vibrational flows was inversely to R and demonstrated a non-linear increase with AVC amplitude (A). Integral viscous dissipation energy for the whole melt volume, which can be processed by integrating Equation (9) over the volume, demonstrated the cubic dependence vs. the oscillating velocity and inverse proportionality to edge radius: P_w ~ A³ω³R⁻¹. The maximum viscous dissipation (Figure 9) near the disc edge is also proportional to oscillating velocity and inverse proportionality to edge radius A³ω³R⁻¹, but it increases faster with 3.8 number exponents.

In the case of NaNO₃ melt, the most stable values of viscous dissipation energy were obtained, since at f = 25 Hz, A = 500 μm its flow, despite the formation of separation vortices, is still stable. Under these conditions, CdTe melt is in a stronger turbulent regime, which is expressed in variable values of the maximum dissipation energy. In the case of Ge-melt at f = 25 Hz, A = 500 μm, an extremely unstable turbulent regime of oscillatory flow is achieved, so we calculated the dissipation energy for the parameters f = 25 Hz, A = 100 μm.

Recalculation of an instant maximum dissipation rate for 0.1 mm edge radius, results in ~220 Joule/mole per second in the case of NaNO₃ melt; and ~330 Joule/mole per second for CdTe melt in one cubic nanometer near the disk edge. These values match the approximation for energy of cluster formation in a liquid phase [21]. For Ge at R = 0.1 mm we are able to introduce only ~2–3 Joule/mole per second in one cubic nanometer near the disc edge. This energy is enough to organize the laminar flow motion in the melt for desired heat-mass transfer, but it is insufficient to change the Ge melt thermodynamic state with respect to components’ composition.

4. Discussion

The resulting dependence of the velocity v on the intensity I (Equation (9)) is primarily a mathematical correlation, but it is also possible to speculate on its physical interpretation. Using the data approximation in Figure 6, we can conduct the linear regression equation in the form:

$$v=b\mathit{I}+v_0=b·\frac{{\rho }^{0.5}{\mathit{A}}^3{\mathit{\omega }}^2}{\mu \mathit{R}}+v_0,$$

(12)

Equation (12) could be transformed to the form:

$$\rho {\left(v-v_0\right)}^2=b^2·\frac{{\mathit{A}}^6{\mathit{\omega }}^4}{{\mathit{R}}^2{\nu }^2}={\left(\frac{{\mathit{A}}^2\mathit{\omega }}{\nu }·\frac{\mathit{A}}{\mathit{R}}·b\mathit{\omega }\right)}^2,$$

(13)

where b is the proportionality coefficient, $\nu =\mu /\rho$ is the kinematic viscosity. We see that the left side of the equation is reduced to the form of dynamic pressure of the melt. On the right side of the equation, the factor ${\mathrm{Re}}_s={\mathit{A}}^2\mathit{\omega }/\nu$ represents the Reynolds number for a stationary melt flow during oscillatory motion of the disk, and the factor $\epsilon =\mathit{A}/\mathit{R}$ is the ratio of the oscillation amplitude to the characteristic size of the body. Also, the factor $\epsilon$ can be interpreted as the ratio of the oscillatory velocity to the flow velocity, and this hints to us that this ratio must contain some other characteristic size, which is a function of R, since the flow velocity near the edge of the disk is 10 times greater than the oscillatory flow velocity, as stated earlier. In the case of direct proportionality, the missing parameter is included in b.

Both quantities are dimensionless, therefore, in the expression:

$$\rho {\left(v-v_0\right)}^2={\mathrm{Re}}_s^2·{\epsilon }^2·{\left(b\mathit{\omega }\right)}^2,$$

(14)

the factor ${\left(b\mathit{\omega }\right)}^2$, obviously, must also be equal in dimension to pressure. Then, you can evaluate:

$$\left[b^2\right]=\left[m/l\right]=\left[\mathrm{k}\mathrm{g}·{\mathrm{m}}^{-1}\right],$$

(15)

and get that:

$$\rho {\left(v-v_0\right)}^2={\mathrm{Re}}_s^2·{\epsilon }^2·\frac{m{\mathit{\omega }}^2}{l}={\mathrm{Re}}_s^2·{\epsilon }^2·∆p.$$

(16)

The resulting expression is very similar in appearance to Bernoulli’s law for an ideal fluid with corrections for viscous dissipation. From the obtained dependence, we can assume that $v_0$ is the limiting speed at which the speed of the secondary vortex along the axis is determined only by the speed of oscillation of the disk. At low frequency and amplitude of oscillations, ${\mathrm{Re}}_s^2$ will decrease quite quickly, and then the right side will tend to zero. Also, an increase in the radius of the disk edge will lead to a decrease in the parameter ${\epsilon }^2$, which in the limiting case will also lead to zeroing of the right side. At small values of ${\mathrm{Re}}_s^2$ and ${\epsilon }^2$, a regime of complete flow around the disk should occur when it oscillates in a liquid, where the maximum velocity is concentrated in a thin layer near the disk. However, to confirm this hypothesis, it is necessary to conduct a study of fluid flow in the region of small amplitudes and vibration frequencies, as well as expand the nomenclature of materials.

To correctly identify the value of $b^2$ and determine what pressure gradient ∆p is implied in the resulting recording, it is necessary to conduct additional studies of the velocity values for sodium nitrate in the laminar region, as well as for germanium in the turbulent region, and to consider the relationship of parameters with other velocities. For instance, the maximum speed at the edge of the disk, as well as in the space between the disk and the crucible. These questions are of interest for further development of the study. Now, we can only assume that the pressure gradient, judging by the analysis of the b dimension, is somehow related to the inertial forces arising when the disk oscillates, as indicated by the fact that we presented ∆p above as:

$$∆p={\left(b\mathit{\omega }\right)}^2=\frac{m{\mathit{\omega }}^2}{l}=\frac{m·l{\mathit{\omega }}^2}{l^2}=\frac{m}{l^2}·\left|\frac{d\overline{v}}{\mathit{dt}}\right|=\frac{\left|\overrightarrow{F}\right|}{S}.$$

(17)

The obtained proportional dependence is useful in that it allows estimating the order of the speed of movement of the melt in the region near the crystallization front, since the flow motion speed along the axis and along the bottom wall are approximately equal for laminar stationary oscillatory motion. From the correlation (Figure 6) we obtain that the slope of the curve is 8.3·10⁻⁵ kg^0.5 m^−0.5. This value, together with the oscillation intensity I, can be used when estimating the oscillatory contribution of melt motion in real crystal growth processes, where the main driving force of the process is thermal convection or diffusion.

In conclusion, it should be noted that the AVC method is very energy efficient in terms of controlling flows in the liquid phase [22,23,24]. At the same time, the AVC method allows one to potentially overcome the scaling problem when increasing the diameter of the growth system from 20 mm to 200 mm or more without significantly changing the design of the vibration actuator.

5. Conclusions

The effects of AVC can be divided into two parts: (1) it is possible to control the oscillatory movement of the liquid in the melt and create vibration flows controlled in speed and configuration; (2) creation of average fast stable flows that provide a given heat and mass transfer and dominate over thermoconvective flows.

We were able to determine the linear dependence of the maximum melt velocity along the central axis on such system parameters as amplitude, frequency, radius of the disk edge, density, and viscosity of the melt when the disk oscillates in a closed crucible. An assessment was made of the intensity of oscillations I, at which a stable mode of oscillatory flow is maintained, which ranges from 10 to 50 kg^−0.5 m^1.5 s⁻¹. An analysis of the obtained dependence was also carried out and a possible relationship between the parameters was identified, however, to confirm the obtained result it is necessary to collect a large sample of data.

When considering thermoconvective flow, for example, in the VGF method, without taking into account the crystallization heat of the substance, we can assume that there will be practically no flow motion in the melt due to the fact that the hotter melt layer will be located above the colder one. In this case, the side heater creates an axial temperature gradient, which leads to the melt movement from bottom to top along the wall. Taking into account the crystallization heat of the substance, the melt will be set in motion. Depending on the crystal growth rate, the crystallization heat value, and the thermal conductivity of a given substance in the melt and solid, thermal convection in the absence of vibrations can be easily assessed [24].

The general analysis we performed in this work let us provide an answer to what vibration parameters must be set in order to form a vibrational melt flow capable of leveling the crystallization front for substances with different physico-chemical and thermodynamic parameters.

Further research involves conducting physical experiments in an isothermal version in order to validate the data obtained. Next, the performance of the model will be tested on non-isothermal systems, taking into account the crystallization process on low-temperature model substances. After this, the developed model can be used to optimize the growth processes of high-temperature crystals.