The Mathematical Model for the Secondary Breakup of Dropping Liquid

: Investigating characteristics for the secondary breakup of dropping liquid is a fundamental scientiﬁc and practical problem in multiphase ﬂow. For its solving, it is necessary to consider the features of both the main hydrodynamic and secondary processes during spray granulation and vibration separation of heterogeneous systems. A signiﬁcant di ﬃ culty in modeling the secondary breakup process is that in most technological processes, the breakup of droplets and bubbles occurs through the simultaneous action of several dispersion mechanisms. In this case, the existing mathematical models based on criterion equations do not allow establishing the change over time of the process’s main characteristics. Therefore, the present article aims to solve an urgent scientiﬁc and practical problem of studying the nonstationary process of the secondary breakup of liquid droplets under the condition of the vibrational impact of oscillatory elements. Methods of mathematical modeling were used to achieve this goal. This modeling allows obtaining analytical expressions to describe the breakup characteristics. As a result of modeling, the droplet size’s critical value was evaluated depending on the oscillation frequency. Additionally, the analytical expression for the critical frequency was obtained. The proposed methodology was derived for a range of droplet diameters of 1.6–2.6 mm. The critical value of the diameter for unstable droplets was also determined, and the dependence for breakup time was established. Notably, for the critical diameter in a range of 1.90–2.05 mm, the breakup time was about 0.017 s. The reliability of the proposed methodology was conﬁrmed experimentally by the dependencies between the Ohnesorge and Reynolds numbers for di ﬀ erent prilling process modes.


Introduction
Scientific research of the processes of vibrational prilling [1,2], granulation [3,4] and separation of gas-dispersed systems [5,6] is an essential problem that has not been wholly studied is the mechanism of the breakup of dropping liquid. In this regard, a crucial hydrodynamic criterion that determines the behavior of droplets in a heterogeneous environment [7] is the Weber number (We) [8]. Its critical value [9] was experimentally evaluated in the articles [10,11], which aimed to define the single droplet breakup and simulations of droplet deformations in airflow.
The main areas of use of dispersed systems are as follows: production of granulated fertilizers; development of compact cooling towers of power plants and monodisperse nozzles; granulation of nuclear fuel; granulation of vitamin preparations; development of express systems for the diagnosis of cells and bioactive substances; development of new composite-based granular materials; designing micro dispensers for medical and biological products; dispensers for rare substances; droplet generators for studying combustion processes, as well as for heat and mass transfer; drip radio space heat exchangers and contactless fueling systems. Consequently, the dispersed systems have gained popularity since they ensure resource conservation, environmental friendliness and quality of new products obtained in technological processes.
When the critical Weber number is exceeded, the secondary breakup occurs. However, the proposed values differ significantly [12]. Mainly, D. Pazhi and V. Galustov in 1984 developed the fundamentals of spraying liquids [13]. As a result, it was experimentally obtained that the Weber number is in a range of We cr = 4-20.
Moreover, under conditions close to critical, the mechanisms of droplet breakup significantly differ. In this case, there are two types of droplet breakup: vibration mode and consequent destruction of the droplet with the formation of a thin film. Additionally, it was experimentally established that the mechanism of breaking up the droplet liquid depends on the hydrodynamic characteristics of the nonstationary flow, which affects the duration of the effect of the gas flow on the droplet.
According to the various conditions of the process, which were previously studied by S. Ponikarov [14] from the Department of Machines and Apparatus for Chemical Production at Kazan National Research Technological University, L. Ivlev and Y. Dovgalyuk [15] from the Research Institute of Chemistry at Saint Petersburgh University, as well as A. Cherdantsev [16] from the Nonlinear Wave Processes Laboratory at Novosibirsk State University, there are different mechanisms of the droplet breakup. The first one is blowing up the middle of a droplet with the subsequent breaking up the toroidal particle. The second one is in the disordered breakup of the droplet into several particles. The last one is in tearing small droplets from the surface of the droplet blown by the flow.
Thus, the study of the characteristics for the secondary breakup of the dropping liquid is an urgent scientific and practical problem. Its solution will allow considering the peculiarities of the operating processes of the prilling [17,18], separation of heterogeneous systems [19,20], pneumatic classification [21,22], spraying of liquid mixtures [23,24] and other hydromechanical, heat and mass transfer processes.
In major technological processes, the breakup of droplets and bubbles occurs with the simultaneous action of several dispersing mechanisms. However, there are significant difficulties in creating mathematical models of breakup liquid droplets and finding accurate analytical solutions. Two types of instability of the typical wave nature occurring in different particle surface parts have been established.
The first mechanism, the Kelvin-Helmholtz instability, occurs in the presence of a shift between the layers of a continuous environment or when two contact media have a significant difference in their velocities [25]. In this case, the boundary layer is failed, when the Weber number is relatively small. When the Weber number exceeds its critical value, microparticles are separated due to wave perturbations on the droplet's side surface.
The Rayleigh instability [26] is the second mechanism related to the spontaneous increase in pressure, density and velocity pulsations in the inhomogeneous environment in a gravitational field or moves with acceleration [27]. For example, on the frontal surface of a falling droplet, oscillations Energies 2020, 13, 6078 3 of 17 occur under the free-fall acceleration. The side surface layer oscillates as a result of the maximum flow velocity. Simultaneously, the particle breakup mechanism under turbulent pulsations' influence has a different character and acts on the particle surface [28].
Thus, droplets' breakup is quite complicated and determined by the ratio of inertia forces, surface tension, viscosity and other factors.

Literature Review
In the article [9], it is noted that a complete model for the aerodynamic breakup of liquid droplets has not been developed. All existing ideas about the laws of breakup and the determining parameters are mainly obtained from the experiment. The most studied is droplet breakup in shock waves. Notably, the main process parameters are the dimensionless criteria of Weber (We), Laplace (La) or Ohnesorge (Oh), Bond (Bo) and Strouchal (Sh). Notably, the Weber number has the most significant influence on the breakup mode [29]. Moreover, according to the Weber number, there are different classifications of breakup modes, distinguishing six [30] or five [9] main mechanisms. Furthermore, according to the previous experience, the most consistent with the up-to-date concepts, the classification of droplet breakup modes for low-viscosity liquids should be considered.
In 1988, M. Clark conceptually studied the droplet breakup model in a turbulent flow [31]. In 1998, W. Ye, W. Zhang and G. Chen [32] investigated the effect of Rayleigh-Plateau instability for a wide range of wavelengths numerically. In 1999, L. Ivlev and Y. Dovgalyuk [15] studied blowing up the middle of a droplet to save its toroidal shape with the subsequent disintegration. They developed the methodology to research the disordered destruction of the droplet. In 2020, G. Chiandussi, G. Bugeda and E. Onate [33] proposed variable shape definition with C 0 , C 1 and C 2 continuity functions in shape optimization problems. In 2019, A. Cherdantsev, in his D.Sc. thesis, "The wave structure of a liquid film and the processes of dispersed phase exchange in a dispersed-annular gas-liquid flow" [16] investigated the process of tearing small particles from the surface of a droplet.
Since the Khyentse-Kolmogorov's equation cannot be applied to describe the breakup of liquid due to the velocity gradient, which excludes the hypothesis of isotropic turbulence. Consequent research of V. Sklabinskyi and B. Kholin solved this problem [34] by considering the velocity gradient.
The breakup's nature is significantly different for different characteristic velocities of the continuous and dispersed phases' relative motion. In 1984, S. Ponikarov, in his D.Sc. thesis, "Droplet breakup in centrifugal equipment of chemical plants" [14] conducted a comparative analysis of theoretical and experimental studies of the droplet breakup process. It was established that there are several fundamental mechanisms of the breakup, which correspond to different ranges of the Weber number.
The viscous friction can crush liquid droplets and bubbles entering the shear flow of a continuous environment. Thus, G. Kelbaliev and Z. Ibragimov, in their research work [35], studied the droplet breakup process in the Quetta stream. They found that a droplet breakup in a turbulent gas flow occurs differently if the environment's density is insignificant compared to the droplet density. In this case, inertial effects play an essential role in the mechanism of the droplet breakup. Experimental confirmation of a critical Weber number was presented in the research work [36].
It has also been established that there are at least two mechanisms of the breakup [37]. At a specific ratio of the length of the drop to its diameter, the breakup occurs with the formation of two new particles of approximately equal size. Suppose this ratio is not met, the droplet thin in several places at once. Another mechanism of the breakup is that the smaller droplet is moved away from a larger one. This mechanism is observed when the vortex velocity becomes critical, making the unstable droplets.
Additionally, the developed mathematical model describes gas-dynamic processes during mixture formation and evaporation of liquid droplets in the nonstationary supersonic flow. This technique allows for designing air-jet engines, power plants, high-performance ejectors, heaters and various technological devices. Additionally, droplet breakup's experimental technique in a stream with a shock wave was presented in [38]. Moreover, the research [39] is devoted to identifying zones of jet self-oscillations as the Hartmann effect and determining the corresponding oscillation frequency at which the process of a fluid jet breakup occurs.
In 2019, Z. Wu, B. Lv and T. Cao improved the Taylor analogy breakup (TAB) and Clark models for droplet deformations [40]. However, these models still operate by linear differential equations, in which inertia, stiffness and damping coefficients are expressed through the similarity criteria. Such approaches, mainly, do not allow substantiating a critical value of the Weber number.
Preliminary analysis of the research works on the process of the breakup of liquid droplets in a wide range of changes in Reynolds, Weber and Laplace numbers [41] indicates the existence of different types of the gas-dynamic breakup differing in intensity and trajectories of detached parts of droplets [42].
Finally, in 2001, B. Gelfand et al. experimentally determined that the Weber number's values less than its critical value do not lead to droplets' breakup [43]. This hypothesis requires a theoretical justification, which is realized below. Moreover, there is still no reliable mathematical model of the breakup process that analytically determines the Weber number's critical value.
The presence of many up-to-date studies also emphasizes the urgency of the considered problem. Particularly, D. Kim and P. Moin [44] developed the subgrid-scale capillary breakup model for liquid jet atomization. X. Li et al. [45] studied the breakup dynamics of the low-density gas-liquid interface during the Taylor bubble formation. A. Salari et al. [46] investigated the breakup of bubbles and droplets in microfluidics. N. Speirs et al. [47] studied jet breakup in normal liquids. C. Tirell, M.-C. Renoult and C. Dumouchel [48] proposed the methodology for measuring extensional properties during a jet breakup. A. Dziedzic et al. [49] studied the substrate effects on the breakup of liquid filaments. J. Zhang et al. [50] investigated the process of an in-fiber breakup. Finally, J.-P. Guo et al. [51] proposed the instability breakup model of fuel jets.
The critical overview of the above studies allows stating the following research gaps. First, since generally, the critical Weber number varies significantly in a range of 4-20 for the secondary breakup of droplets, its range should be narrowed for the particular range of droplet diameters. Second, the proposed methodologies generally operate with similarity criteria and empiric constants. As a result, the droplet breakup time and the Weber number's critical values were not still substantiated analytically.

A Mathematical Model
At the foundation of the breakup of liquid droplets is a mechanism according to which the deformation of a droplet takes the form of an elongated ellipsoid with the subsequent transformation of its shape and disintegrating into two approximately equal particles. Therefore, the process of breakup of a droplet into two equal parts of the spherical shape is considered below.
In the initial stage, this process can be described as a deviation from the equilibrium state. The last one is determined by the action of inertia forces m 1 a C1 , gravity G, Archimedes force F A and surface tension force F σ (Figure 1).
For the description of the proposed mathematical model, the following changes were considered in the geometry of a droplet and the description of the forces' impact on it. At the initial time (t = 0), the drop is spherical. Before its crushing, forces of inertia, gravity and Archimedes act. Further, at the initial stage of secondary crushing, the surface tension force is changed due to the change in the angle of inclination of the tangent to the forming secondary droplet. This angle changes as the drop are breaking.
Consequently, the surface tension force makes a time-varying contribution to the overall force action. In this case, the total mass of the main drop and the satellite remains unchanged. Finally, the secondary breakup time is determined from the separation of a single droplet into a couple of droplets. Under these forces' actions, the mass centers of a couple of mass m1 (kg) reached acceleration aC1 (m/s 2 ). In this case, the fundamental equation describing the motion of a single droplet has the following form: where m1 = ρp(V / − Vs) is the mass of a droplet, kg; aC1 = 0.5d 2 Δz/dt 2 is the acceleration of the mass center (m/s 2 ), determined as the second derivative from the relative distance of the droplet parts Δz (m) with respect to time t (s); Δα is wetting angle (rad).
The following dependencies determine the acting forces mentioned above: where g-gravitational acceleration (m/s 2 ); mp-mass (kg) of a droplet with volume V (m 3 ); ρp, ρm-the densities of the dropping liquid and the environment, respectively (kg/m 3 ); σ-surface tension coefficient (N/m); L-the perimeter of the wetted contour (m). Droplet volume V = 2(V / − Vs) is defined as the double difference of its spherical part (V / ) and segment (Vs): where r = (R − Δr)-the current radius of droplets (m), which is defined as the difference between its initial value R and the magnitude of the change in radius Δr due to its breakup; h = (r − 0.5Δz)-the height of the spherical segment, which is associated with the relative displacement Δz of the parts. The dependence between the change in radius Δr of a droplet and its displacement Δz can be determined from the conservation law for mass: where V0 = 4πR 3 /3-initial volume (m 3 ). Considering Equation (3) after identical transformations allows writing cubic equation concerning dimensionless ratio Δz/R: Under these forces' actions, the mass centers of a couple of mass m 1 (kg) reached acceleration a C1 (m/s 2 ). In this case, the fundamental equation describing the motion of a single droplet has the following form: where m 1 = ρ p (V / − Vs) is the mass of a droplet, kg; a C1 = 0.5d 2 ∆z/dt 2 is the acceleration of the mass center (m/s 2 ), determined as the second derivative from the relative distance of the droplet parts ∆z (m) with respect to time t (s); ∆α is wetting angle (rad).
The following dependencies determine the acting forces mentioned above: where g-gravitational acceleration (m/s 2 ); m p -mass (kg) of a droplet with volume V (m 3 ); ρ p , ρ m -the densities of the dropping liquid and the environment, respectively (kg/m 3 ); σ-surface tension coefficient (N/m); L-the perimeter of the wetted contour (m). Droplet volume V = 2(V / − Vs) is defined as the double difference of its spherical part (V / ) and segment (Vs): where r = (R − ∆r)-the current radius of droplets (m), which is defined as the difference between its initial value R and the magnitude of the change in radius ∆r due to its breakup; h = (r − 0.5∆z)-the height of the spherical segment, which is associated with the relative displacement ∆z of the parts. The dependence between the change in radius ∆r of a droplet and its displacement ∆z can be determined from the conservation law for mass: where V 0 = 4πR 3 /3-initial volume (m 3 ). Considering Equation (3) after identical transformations allows writing cubic equation concerning dimensionless ratio ∆z/R: Energies 2020, 13, 6078 6 of 17 At the initial stage of the drop breakup, when the relative displacement is insignificant in comparison to the initial radius (∆z << R), the above equation is reduced to a linear one, and the dependence between the parameters ∆z and ∆r takes the following form: From the condition of the breakup process completion, it can be established that the reduction of the radius ∆r does not exceed its maximum value ∆r max = (1 − 2 1/3 ) ≈ 0.2. In this case, the last expression can be reduced to a simplified one in a linearized form as ∆z ≈ 0.25∆r ( Figure 2).
Energies 2020, 13, x FOR PEER REVIEW 6 of 17 At the initial stage of the drop breakup, when the relative displacement is insignificant in comparison to the initial radius (Δz << R), the above equation is reduced to a linear one, and the dependence between the parameters Δz and Δr takes the following form: From the condition of the breakup process completion, it can be established that the reduction of the radius Δr does not exceed its maximum value Δrmax = (1 − 2 1/3 ) ≈ 0.2. In this case, the last expression can be reduced to a simplified one in a linearized form as Δz ≈ 0.25Δr ( Figure 2).  (2) dependencies between the dimensionless reduction of the droplet radius (Δr/R) and the dimensionless ratio of the relative displacement of its parts (Δz/R).
The wetting angle Δα included in the droplet motion's fundamental equation is determined by the following trigonometric expression: It is simplified, which occurs at the initial stage of the droplet breakup (considering the relatively small values of parameters Δz and Δr) and applying the first limit.
Additionally, in the first approximation, the perimeter of the wetting contour can be determined by the following dependence: Substitution of expression (2) to Equation (1) considering Formulas (7) and (8) allows rewriting the equation of motion in the following form: Notably, in stationary mode (d 2 Δz/dt 2 = 0; Δz = Δz0 = const), the droplet breakup occurs under the following condition: In this regard, limit value dcr = 2R of the unstable droplets' diameter is determined from the condition of zero displacement Δz0. Thus, in the case of relatively heavy drops (ρp >> ρm), it can be obtained as follows:  (2) dependencies between the dimensionless reduction of the droplet radius (∆r/R) and the dimensionless ratio of the relative displacement of its parts (∆z/R).
The wetting angle ∆α included in the droplet motion's fundamental equation is determined by the following trigonometric expression: It is simplified, which occurs at the initial stage of the droplet breakup (considering the relatively small values of parameters ∆z and ∆r) and applying the first limit.
Additionally, in the first approximation, the perimeter of the wetting contour can be determined by the following dependence: Substitution of expression (2) to Equation (1) considering Formulas (7) and (8) allows rewriting the equation of motion in the following form: Notably, in stationary mode (d 2 ∆z/dt 2 = 0; ∆z = ∆z 0 = const), the droplet breakup occurs under the following condition: Energies 2020, 13, 6078 7 of 17 In this regard, limit value d cr = 2R of the unstable droplets' diameter is determined from the condition of zero displacement ∆z 0 . Thus, in the case of relatively heavy drops (ρ p >> ρ m ), it can be obtained as follows:

A particular Solution of the Model
Applying the method of variations allows considering the nonstationary relative motion of a droplet near its equilibrium position. In this case, the displacement of the droplet ∆z = (∆z 0 + δz) is defined as a sum of the stationary component ∆z 0 and variation δz relative to the equilibrium state, and the differential equation of motion in variations takes the following form: A general solution has the following form: where the following oscillatory parameter is introduced (s −1 ): Integration constants C 1 , C 2 are determined from the initial conditions. Particularly, for zero initial deviation δz(0) = 0 and velocity (dδz/dt)| t = 0 = v 0 it can be obtained C 1 = v 0 /λ and C 2 = 0. Thus, solution (13) takes the following form:

A Critical Value of the Weber Number
The resulting dependence allows evaluating the time T s (s) of droplet breakup. According to the condition: the following time can be obtained: where We = ρ p v 0 2 R/σ is the Weber number as the ratio of the specific inertia and the surface tension forces.
It should be noted that the critical value of the Weber number We cr = 3π/4 1/3 ≈ 5.94 obtained from this formula allows considering individual cases of solving the Equation (17). Notably, in the case of significant velocities (We >> We cr ), the breakup time does not depend on the Weber number due to the relative smallness of the surface tension forces compared with the inertia forces. In this case, the last formula is approximately equal to the following one: Energies 2020, 13, 6078 8 of 17 Introduction of the Ohnesorge and Laplace numbers [52]: allows obtaining the critical values of the Ohnesorge and Laplace numbers: The theoretical line of this dependence is graphically presented in Figure 3.
Energies 2020, 13, x FOR PEER REVIEW 8 of 17 allows obtaining the critical values of the Ohnesorge and Laplace numbers: The theoretical line of this dependence is graphically presented in Figure 3. On the other hand, for relatively slow droplets (We << Wecr), the following approximation can be written:

Characteristics of the Secondary Breakup
In the case of deceleration of a droplet (d 2 Δz/dt 2 < 0), the vibration frequency, at which the resonance of the droplet liquid occurs with its subsequent breakup, is equal to the following one (Hz): This dependence allows obtaining the expression for the critical diameter of droplets under the external vibrational impact: Thus, an increase in vibration frequency leads to a decrease in the size of droplets.

Consideration of the Oscillating Wall
It should be noted that the above-mentioned equation of motion for the liquid droplets does not consider the vibration force Fv (N), which is determined by the dependence previously obtained by L. Blekhman   On the other hand, for relatively slow droplets (We << We cr ), the following approximation can be written:

Characteristics of the Secondary Breakup
In the case of deceleration of a droplet (d 2 ∆z/dt 2 < 0), the vibration frequency, at which the resonance of the droplet liquid occurs with its subsequent breakup, is equal to the following one (Hz): This dependence allows obtaining the expression for the critical diameter of droplets under the external vibrational impact: Thus, an increase in vibration frequency leads to a decrease in the size of droplets.

Consideration of the Oscillating Wall
It should be noted that the above-mentioned equation of motion for the liquid droplets does not consider the vibration force F v (N), which is determined by the dependence previously obtained by where v 0 -amplitude of vibration velocity (m/s); z 0 -the distance from a droplet to the source of vibrations (m). Notably, the vibrating force acts on the particles phase in a continuous medium near an impenetrable wall in the normal direction. The effect of vibrational weighing is caused by the pressure averaged over the oscillation period [54] due to convective inertia forces.
The introduction of this force allows representing the differential Equation (9) in the more generalized form: ∆ ..
In variations δz relative to the stationary position z 0 , the last equation takes the following form: Thus, the droplet breakup condition under the vibrational impact (F v > F σ ) at the stability mode limit corresponds to a zero value of the coefficient before variation δz. This fact allows determining the critical size of droplets under the impact of the vibrational force: This expression indicates that even highly dispersed particles can be involved in the secondary breakup process when the particles approach the oscillating wall.

Consideration of the Resistance Force
For considering the impact of the resistance force F r (N) on the secondary breakup, the right part of the fundamental Equation (9) should be supplemented by the appropriate component. Notably, in the case of conditionally linear dependence of the resistance force on the droplet velocity [55]: the equation of motion in variations takes the following form: where µ m -dynamic viscosity of the environment (Pa·s); ν m = µ m /ρ m -kinematic viscosity of the environment (m 2 /s). The last differential equation has a general solution: where p 1,2 = −n ± √ λ 2 + n 2 are the characteristic equation's roots containing the damping factor n = 4.5πν m /R 2 (s −1 ).
For the above initial conditions, the integration constants are determined as follows: Finally, solution (30) takes the following form: Notably, the previously obtained solution (15) is a particular case (for n = 0) of a more general solution (30).
Considering the resistance force, the time T s of droplet breakup reduced to finding the positive root of the transcendental equation: which is for the case of n << λ with sufficient for practical needs accuracy can be simplified to the following cubic equation: Its positive root can also be determined using the methods for parameter identification of hydromechanical processes [56], particularly by the following iteration procedure: where T s1 -the breakup time determined by the formula (18); i-a number of the current iteration. It is convenient to choose the initial value of the breakup time at the initial iteration as T s <0> = T s1 .

Experimental Research Data
Investigation of hydrodynamic parameters for the secondary droplet breakup was carried out on the modernized vibrational priller of the Department of Chemical Engineering at Sumy State University for the following parameters: the wall thickness-1 mm; hole diameter-1 mm; modeling liquid-water at temperature 20 • C; dispersion medium-air at temperature 25 • C. The design scheme of the experimental stand is given in Figure 4

Experimental Research Data
Investigation of hydrodynamic parameters for the secondary droplet breakup was carried out on the modernized vibrational priller of the Department of Chemical Engineering at Sumy State University for the following parameters: the wall thickness-1 mm; hole diameter-1 mm; modeling liquid-water at temperature 20 °C; dispersion medium-air at temperature 25 °C. The design scheme of the experimental stand is given in Figure 4. The stand functions as follows: Due to the created vacuum by pump 10, the liquid from the bottom of the buffer tank 9 is transferred into the pipeline. After, it is fed into pipe 3 of the upper part of the vibrating priller. Fluid flow is monitored by an electromagnetic flowmeter 12 and regulated by The stand functions as follows: Due to the created vacuum by pump 10, the liquid from the bottom of the buffer tank 9 is transferred into the pipeline. After, it is fed into pipe 3 of the upper part of the vibrating priller. Fluid flow is monitored by an electromagnetic flowmeter 12 and regulated by a valve 11. Through pipe 3, the liquid is fed into the annular collector 4. After this, the operating fluid passes through the filter element 5 as a perforated cylinder. As it passes through it, air bubbles are released from the liquid volume. The filter grid is fixed on the perforated cylinder to prevent the clogging of holes in a basket 1. The liquid then flows to the perforated bottom 1, gradually filling its volume. The filling level is monitored by a liquid level meter 20. Under pressure created by the hydrostatic liquid level, the liquid flows out of holes in the perforated bottom 1.
Using a particular program, PC 13 generates a vibrational signal. Using the low-frequency amplifier 14, the magnetic actuator 6 through the rod 7 brings the resonator 8 into oscillating motion with a given frequency (Hz) and amplitude (V). The resonator disk is placed above the central part of the perforated bottom. The disk and the bottom's liquid gap provide a hydrodynamic connection between all the studied hydromechanical system elements. When imposing the vibrations, the resonator disk performs reciprocating motion. Oscillatory waves propagate as elastic deformations in the liquid and are transferred to the perforated bottom. As a result, regular perturbations are superimposed on the fluid flowing from the holes. This effect causes the breakup of a liquid jet into droplets in places of narrowing.
The obtained images (Figure 5a) can be processed using both the "Matlab" [57] and the "Digimizer" software, particularly by the algorithm for the detection of circular elements [58]. This approach allows determining droplets' size and distances between them (Figure 5b).
in the liquid and are transferred to the perforated bottom. As a result, regular perturbations are superimposed on the fluid flowing from the holes. This effect causes the breakup of a liquid jet into droplets in places of narrowing.
The obtained images (Figure 5a) can be processed using both the "Matlab" [57] and the "Digimizer" software, particularly by the algorithm for the detection of circular elements [58]. This approach allows determining droplets' size and distances between them (Figure 5b). It should be noted that the camera was positioned perpendicular to the screen at a distance of 500 mm. This is the minimum focal length from the subject to the focal plane during photography. The photographs of the obtained drops were processed using the program Digimizer that automatically detects an object by measuring the geometrical characteristics. A marking grid is required to set the measurement scale (pixels/mm) for the program when processing photos. Measurement results were obtained using built-in algorithms to detect objects, particularly for finding round shapes and sizing them.
The experimental results for different liquid levels are summarized in Table 1 and graphically presented in Figure 6.  It should be noted that the camera was positioned perpendicular to the screen at a distance of 500 mm. This is the minimum focal length from the subject to the focal plane during photography. The photographs of the obtained drops were processed using the program Digimizer that automatically detects an object by measuring the geometrical characteristics. A marking grid is required to set the measurement scale (pixels/mm) for the program when processing photos. Measurement results were obtained using built-in algorithms to detect objects, particularly for finding round shapes and sizing them.
The experimental results for different liquid levels are summarized in Table 1 and graphically presented in Figure 6.  The data in Table 1 presents the mean values. Notably, the standard deviation does not exceed the value of 0.02 mm.

An Example of Numerical Calculation
Based on the proposed methodology, an example of numerical calculation can be realized. Particularly, as initial data, the following parameters were chosen: the modeling liquid was water at a temperature of 20 °C; density of particles ρp = 1·10 3 kg/m 3 ; surface tension coefficient σ = 7.3·10 -2 N/m; kinematic viscosity of the medium νm = 1.5·10 -6 m 2 /s; gravitational acceleration g = 9.81 m/s 2 . The oscillation frequency of the actuator f = 260 Hz, at which the secondary breakup occurs with the formation of satellites of approximately the same size compared to the diameter of the main droplet; the critical value of the Weber number Wecr = 5.94; the experimentally obtained range of droplet diameters is in a range of 1.90-2.05 mm ( Table 1).
The values of critical diameters determined by Formulas (11) and (26) exceed the holes' diameter, consequently were not considered. Thus, for this experimental case, the expression (22) allows evaluating a breakup droplet diameter of 1.98 mm. This value falls within the range of observed values of droplet diameters with a relative error of about 4%.

Discussion
Due to the analytical and experimental results described above, the following detailed discussion should emphasize a range of applicability of the proposed mathematical model and substantiate its advantages compared to the well-known models. Remarkably, the Rayleigh instability model describes why and how a falling stream of fluid breaks up into satellites with less surface area suitable for substantiation of polydisperse modes. Additionally, according to the C 0 continuity model [59], the time The data in Table 1 presents the mean values. Notably, the standard deviation does not exceed the value of 0.02 mm.

An Example of Numerical Calculation
Based on the proposed methodology, an example of numerical calculation can be realized. Particularly, as initial data, the following parameters were chosen: the modeling liquid was water at a temperature of 20 • C; density of particles ρ p = 1 × 10 3 kg/m 3 ; surface tension coefficient σ = 7.3 × 10 −2 N/m; kinematic viscosity of the medium ν m = 1.5 × 10 −6 m 2 /s; gravitational acceleration g = 9.81 m/s 2 . The oscillation frequency of the actuator f = 260 Hz, at which the secondary breakup occurs with the formation of satellites of approximately the same size compared to the diameter of the main droplet; the critical value of the Weber number We cr = 5.94; the experimentally obtained range of droplet diameters is in a range of 1.90-2.05 mm ( Table 1).
The values of critical diameters determined by Formulas (11) and (26) exceed the holes' diameter, consequently were not considered. Thus, for this experimental case, the expression (22) allows evaluating a breakup droplet diameter of 1.98 mm. This value falls within the range of observed values of droplet diameters with a relative error of about 4%.

Discussion
Due to the analytical and experimental results described above, the following detailed discussion should emphasize a range of applicability of the proposed mathematical model and substantiate its advantages compared to the well-known models. Remarkably, the Rayleigh instability model describes why and how a falling stream of fluid breaks up into satellites with less surface area suitable for substantiation of polydisperse modes. Additionally, according to the C 0 continuity model [59], the time of droplet breakup is determined more simplistically, like expression (18). However, a more general solution (17) clarifies this approach because the droplet breakup's mathematical model is a nonlinear one.
Additionally, it should be noted that the numerical simulation models commonly use two spray breakup models. The first one is the TAB model [32]. However, it should be pointed out that the TAB model is recommended for only low-Weber-number injections and is well suited for low-speed sprays into a standard atmosphere. Second, for Weber numbers greater than 100, the wave model is more applicable [33]. Consequently, the proposed analytical model allows closing the gap between these two extreme cases. In addition, the proposed method does not operate only with similarity criteria. It is precisely operated by the initial equations reflecting the physical and geometrical essence of the secondary breakup. As a result, the proposed approach's main advantage is the possibility of analytical substantiation of the Weber number's critical value in a range experimentally predicted by previous researchers. Moreover, it should be noted that according to the TAB model, the breakup time is determined numerically using empiric constants, which are dependent on the Ohnesorge, Taylor and Weber numbers. On the other hand, the proposed technique based on Equations (17) and (32) makes it possible to estimate droplet breakup time analytically.
It should also be noted that in the case of mainly used prilling equipment, the atomization mode can be obtained for more than 800 Hz. This frequency depends on the height of a liquid layer, a medium's physical properties and the holes' geometry. In addition, droplet breakup time can be determined using the particle image velocimetry (PIV) method [60].
Moreover, the proposed model allows determining Weber's critical values, Ohnesorge and Laplace numbers for different Reynolds numbers. Particularly, expressions (21) describe the dependencies for the Ohnesorge and Laplace numbers' critical values for the Reynolds number's different values. It should be noted that the obtained range is additionally proved by the experimental results data [61].
Finally, this article significantly supplements the existing mathematical models. It allows predicting the modes of liquid jets breakup and developing new equipment for granulating products with improved characteristics.

Conclusions
Thus, the article investigates the secondary breakup of dropping liquid of the dispersed phase in chemical technology processes. As a result, a mathematical model of the nonstationary decay of droplets was developed. This model considers the impact of volume and surface forces on the relative displacement of the decaying droplet.
The simulation results have allowed developing the dependence for a droplet's critical size and determining the required vibrational frequency. Additionally, the critical diameter of unstable droplets was obtained, and the dependence of breakup time was calculated. These data have allowed analytically determining the critical value of the Weber number.
Additionally, the dependencies between the Ohnesorge, Laplace and Reynolds numbers are obtained analytically and proved experimentally. Particularly, for droplet diameters in a range of 1.6-2.6 mm, the Reynolds and Ohnesorge numbers are in a range of Re = (3.1-7.2) × 10 3 and Oh = (0.35-0.81) × 10 −3 , correspondently.
The reliability of the achieved results was confirmed by the fact that the critical value of the Weber number We cr = 5.9 is in a range We cr = 4-20, obtained experimentally by the previous researchers D. Pazhi and V. Galustov, as well as by consistency between the determined critical size of decaying droplets 1.98 mm and the experimentally obtained range of diameters 1.90-2.05 mm for the droplet breakup at the oscillation frequency 240 Hz for the example of the modernized vibrational priller of Sumy State University. In this case, the relative error is about 4%.
The obtained results can help create appropriate techniques and methodologies for designing the vibrational prillers to get monodispersed prills and granules. Moreover, the proposed mathematical model can also be extended to ensure the gas-dynamic equipment's reliability for the separation of gas-liquid systems.