Numerical Simulation of Acoustic Agglomeration Using the Multi-Monte Carlo Method

1. Introduction

The presence of ambient particulate matter (PM) in the atmosphere has been identified as a serious problem having adverse effects on both human health and the environment [1]. Increased PM concentration in the atmosphere will reduce visibility, but can also affect ecosystems like plants or soil water via its deposition [2].

Based on their size, inhalable PM can be categorized into coarse PM10 (having aerodynamic diameter, d, from about 2.5 μm up to 10 μm), fine PM2.5 (d equal to or less than 2.5 μm) and its ultrafine subset (d less than 0.1 μm). All of these fractions have been found to be toxic by recent epidemiological studies [3]; however, the fine and ultrafine components, which are mainly produced by combustion activities, are even more dangerous [4,5], due to two main factors. First, they consist of carbonaceous combustion particles with chemical components like polycyclic aromatic hydrocarbons (PAHs) that have increased carcinogenic potential [6]; second, because of their size, they can penetrate deeper into the human lungs and cause respiratory (chronic bronchitis, asthma or lung cancer) and cardiovascular diseases (cardiac arrhythmias or heart attacks) [7]. In 2021, the World Health Organization (WHO) revised its Air Quality Guidelines (AQG) with more stringent allowable limits for annual average PM2.5 concentrations compared to the earlier version of 2005 (to 5 μg/m³, down from 10 μg/m³), thus recognizing the severe effect of PM on public health [8,9]. Additionally, in May 2024 the US Environmental Protection Agency (EPA) reduced, for the first time in over eight years, its annual health-based exposure standard for PM2.5 from 12.0 µg/m³ to 9.0 µg/m³, aiming to prevent premature mortality attributed to atmospheric pollution [10]. In an effort to harmonize with the WHO, the European Ambient Air Quality Directive was also revised at the end of 2024 [1]. Even though PM can originate from natural sources (fugitive dust and sand, natural gaseous emissions, airborne pollen, forest fires, sea salt spray) [11], its main sources are anthropogenic emissions such as burning fossil fuels for transportation [12], power generation [13], indoor biomass burning [14], construction sites [15], wildfires and managed burns [16].

Regarding the particulates emitted from the power sector, coal power plants were traditionally the major source of pollutants [17]. Although there is a concerted worldwide effort to shut down coal power plants to achieve the Paris Agreement goals [18], countries like China [19], India [20], South Africa [21], and Indonesia [22] still rely heavily on coal-fired power plants for electricity, despite being major coal producers.

Filtering devices such as electrostatic precipitators, baghouses, wet scrubbers, cyclone dust collectors, and granular bed filters have been employed in power plants for years to capture particles emitted from the combustion of fossil fuels [23]. While these devices exhibit an efficiency exceeding 99.9% for particles with a diameter larger than 2.5 μm [24], their retention capacity diminishes for fine and sub-micron particles [25]. To improve the efficiency of these devices, one approach is the preconditioning of the exhaust gas stream in order to increase the size of the particles with the aid of various agglomeration technologies, like heterogeneous condensation [26] or the application of an external force field [27], such as acoustic [28] or magnetic [29] fields. There is also the chemical agglomeration process, where chemical agents are used to assist the growth of smaller particles that can be easily captured [30]. Using acoustic fields to promote particle agglomeration is the focus of the current work.

Figure 1. Schematic concept of the acoustic agglomeration preconditioning of solid particles found in industrial gas emissions.

Figure 1. Schematic concept of the acoustic agglomeration preconditioning of solid particles found in industrial gas emissions.

During acoustic agglomeration preconditioning, high-intensity acoustic fields are applied to the flue gas stream, inducing interaction effects among suspended particles, leading to collisions and agglomeration [31]. Testing of this methodology within industrial-scale chambers has demonstrated a potential reduction of up to 70% in the number concentration of micron-sized fly ash particles [32].

The observations and studies of acoustically induced particle interactions have a considerable historical precedent: William Ostwald was the first researcher to suggest as early as 1880 that sound waves could be used as a means to collect small droplets; nevertheless, he did not conduct any acoustic-agglomeration-specific experiments. The earliest documented instance of sound waves affecting droplet size is credited to an anonymous reader of a Norwegian newspaper. This individual observed a disruption in the color pattern of the rainbow during a thunderstorm. The observer theorized that the observed anomaly could be attributed to the acoustic phenomenon of thunder, given the near-simultaneous nature of the events [33]. His assumption was supported later by both theoretical explanations [34] and numerical experiments [35,36].

By the conclusion of the 1920s, independent investigations conducted by researchers Wood and Loomis [37] and Andrade and Lewer [38] yielded a consistent observation: particulate matter subjected to acoustic waves demonstrated a propensity for clustering and aggregation. Subsequently, the phenomenon of acoustic agglomeration became the focus of numerous experimental studies, including those by Patterson and Cawood [39], Andrade [40] and St. Clair [41]. Post-World War II, these investigations were continued in Russia and Poland, resulting in significant contributions by Mednikov [42] and Shirokova [43]. Thereafter, numerous studies, both theoretical and experimental, have been performed on the subject of acoustic agglomeration, targeting the retention of PM from coal-fired ash particles [28,31,32,44,45,46,47,48,49,50,51,52,53,54,55].

Many physical mechanisms are believed to contribute to the acoustic agglomeration process: the first one to be accounted for is the orthokinetic mechanism, where collisions between particles of different sizes are attributed to the different ways they respond to the acoustic oscillations of the fluid [56]. The variance in particle motion within the acoustic field, described by each particle's entrainment factor, creates a relative velocity that facilitates inter-particle collisions [57]. The orthokinetic mechanism, while considered the primary driver of acoustic agglomeration, fails to account for collisions between particles of uniform size; a phenomenon substantiated by experimental observations [43,58]. Consequently, the hydrodynamic mechanism has been proposed, positing that inter-particle interactions, even among identically sized particles, are induced by hydrodynamic forces generated during close proximity. Two principal hydrodynamic mechanisms are discernible: the first arises from an effect known as the mutual radiation pressure effect (MRPE) and the second from the acoustic wake effect (AWE). In the instance of mutual radiation pressure, hydrodynamic forces known as Bjerknes forces originate from the nonlinear interactions between the particle scattering wave and the incident wave, recognized as secondary acoustic forces [59,60]. Conversely, concerning the acoustic wake effect, the presence of a particle within the asymmetric flow field surrounding another moving particle leads to a decrease in drag force and an acceleration toward the other particle, culminating in their collision over successive acoustic cycles [61]. In addition to the aforementioned primary mechanisms, the literature suggests that particle interactions may be influenced by additional acoustic phenomena such as the primary acoustic radiation force [62], acoustic streaming [63], and acoustic turbulence [64]. However, these influences are considered small relative to orthokinetic and hydrodynamic interactions and are consequently often disregarded in agglomeration studies [65]. All these mechanisms collectively underlie the dynamics of particle convergence and collision. Subsequent to collision, agglomeration ensues due to prevailing attractive forces, specifically van der Waals forces [66].

Numerical modeling and simulation are valuable tools for investigating the influence of various mechanisms on the overall efficiency of acoustic agglomeration in enhancing particle size. The standard approach to numerically analyzing acoustic agglomeration dynamics involves solving the PBE [67] or its discrete form, the Smoluchowski equation [68]. The solution of this equation depicts the temporal evolution of a particle size distribution, contingent upon the accurate definition of the collision rate. It can be solved through deterministic methodologies, such as the method of moments [69], the sectional method [70], the discrete method [71], the discrete-sectional method [72], or stochastic techniques commonly referred to as Monte Carlo (MC) methods [73]. In contrast to deterministic integration, stochastic methodologies employ MC simulations to model the dynamic behavior of a finite ensemble of particles and track their history over time [67].

Figure 2. Schematic illustrations of the three major acoustic agglomeration mechanisms.

Figure 2. Schematic illustrations of the three major acoustic agglomeration mechanisms.

The objective of the present study is to investigate numerically the acoustic agglomeration of liquid particles. Experimental data on glycol fog droplets agglomerated in an acoustic chamber are used to validate the MMC method. A key contribution of this approach is the incorporation of the acoustic wake effect into the agglomeration model employed for the simulations.

2. Methods and Models

The PBE, accounting solely for agglomeration phenomena, is expressed by the following nonlinear integro-differential equation, originally formulated by Friedlander [74]:

$${\displaystyle \frac{\partial n_p\left(\upsilon ,t\right)}{\partial t}}={\displaystyle \frac{1}{2}}{\int }_0^{\upsilon }\mathsf{\beta }\left(\upsilon -u,u\right)n_p\left(\upsilon -u,t\right)n_p\left(u,t\right)du-{\int }_0^{\mathrm{\infty }}\mathsf{\beta }\left(\upsilon ,u\right)n_p\left(\upsilon ,t\right)n_p\left(u,t\right)du-R\left(\upsilon \right)n_p\left(\upsilon ,t\right)$$

(1)

where n_p(υ,t) is the particle size distribution function at time t and β(υ,u) is the agglomeration frequency function of two particles of volumes υ and u. The first term on the right-hand side of Equation (1) is the birth term, accounting for the formation of a particle of volume υ due to an agglomeration event between a particle of volume u and a particle of volume (υ-u). The factor of 1/2 arises from the fact that each agglomeration event involves two particles, necessitating a correction to avoid double counting. The subsequent term represents the attrition of particles with volume υ due to agglomeration with other particles, known as the death term. The third term accounts for particle loss through deposition, as dictated by the deposition kernel, R(υ). The left-hand side of Equation (1) denotes the temporal variation in the number concentration of particles with volume υ. This formulation assumes a spatially homogeneous system, thereby restricting the analysis to the temporal evolution of the particle size distribution function [75]. The term β(υ,u), also referred to as the collision frequency function or kernel, encapsulates the physical mechanisms that drive inter-particle collisions and agglomeration and quantifies the rate at which two particles of volumes υ and u collide. The following sections will provide an overview of prevalent kernels and detail the numerical methodology employed to solve the PBE.

MC methodologies can be categorized into two distinct classes based on the temporal discretization approach: time-driven and event-driven. In time-driven simulations, a predetermined time increment is specified, and all feasible events occurring within this increment are executed. Conversely, event-driven simulations involve the immediate implementation of an event, followed by an adjustment of the simulation time by a corresponding duration [67]. MC methods can be further classified into two primary classes based on the temporal evolution of the simulation volume or the simulation particle number. The first is known as the “constant volume” approach and the second as the “constant number” approach. During the constant volume simulation, the control volume remains invariant despite fluctuations in the number of particles. However, adjustments are made at discrete time intervals when the particle count surpasses or falls below predefined thresholds, thereby ensuring the statistical independence of the solution [73,76]. On the other hand, during the constant number simulation, the number of particles remains constant. Thus, new particles are continuously added or removed from the initial ensemble, contingent on the dynamic process, to uphold the constant particle number [77].

The Multi-Monte Carlo (MMC) methodology, utilized in the present work for the numerical solution of the PBE, was initially developed to reduce computational cost and improve upon the accuracy of the Direct Simulation Monte Carlo (DSMC) method [78]. In addition to operating under a “constant volume and number” regime, the MMC methodology incorporates a “time-driven” approach for time-step management [79]. This time-driven approach entails tracking a simulation particle over a discrete time interval, considering all potential agglomeration events within this period. This consideration is predicated on the assumption that events are independent when the time-step is sufficiently small [80]. A key distinction from the DSMC methodology lies in the representation of simulation entities. Rather than representing individual particles, these entities represent groups of particles with similar properties, referred to as "fictitious" particles. While the concept of fictitious particles has been previously employed in constant number MC methods, the MMC methodology uniquely features a dynamic weight assignment. This weight, denoting the number of real particles represented by a simulation particle, varies in response to prevailing dynamic processes. Specifically, agglomeration reduces the weight and increases the size of simulation particles, whereas deposition only decreases the weight [81].

Within the MMC methodology, the principal stages deal with the determination of each fictitious particle's weight, the establishment of the time step, the identification of a potential agglomeration partner, the evaluation of an agglomeration event incidence, and the subsequent processing of particle agglomeration. Initially, the adjustment of the weight value for each simulation particle is required. Before this adjustment, it is necessary to ascertain the number of discrete size classes that will accurately represent the overall particle distribution, as well as the total number of simulation particles to be employed. Subsequently, the weight of each simulation particle is computed [82]:

${\left(kwt\right)}_i={\displaystyle \frac{\mathrm{Number}\mathrm{of}\mathrm{real}\mathrm{particles}\mathrm{in}\mathrm{class}i}{\mathrm{Number}\mathrm{of}\mathrm{fictitious}\mathrm{particles}\mathrm{in}\mathrm{class}i}}$ (2)

To ensure accurate tracking of all fictitious particles and potential events through appropriate time step selection, it is necessary to avoid the occurrence of multiple events within a single time step. This criterion guarantees the enumeration of each individual event. This stipulation necessitates referencing the dynamic events of actual particles and utilizing kernel information about their frequency. Specifically, the rate of real agglomeration events between particles of the ith-group and the jth-group, per unit time and per unit volume V, is derivable from the equation provided below [78]:

${\mathsf{\varphi }}_{ij}={\mathsf{\beta }}_{ij}\times {\displaystyle \frac{{\left(kwt\right)}_i}{V}}\times {\displaystyle \frac{{\left(kwt\right)}_j}{V}}$ (3)

The agglomeration rate for a representative particle from the ith-group interacting with all other particles can be determined by:

$C_i={\displaystyle \frac{1}{V^2}}\left[{\sum }_{j=1,j\ne i}^{N_f}\left({\mathsf{\beta }}_{ij}\times {\left(kwt\right)}_j\right)\right]\approx {\displaystyle \frac{1}{V^2}}\left[{\sum }_{j=1}^{N_f}\left({\mathsf{\beta }}_{ij}\times {\left(kwt\right)}_j\right)\right]$ (4)

Consequently, the temporal duration separating two successive agglomeration events of a particle may be expressed as:

$t_{i,coag}^{fict}={\displaystyle \frac{1}{V{C}_i}}$ (5)

and the time step is calculated as:

$\mathsf{\Delta }t\le {\mathrm{m}\mathrm{i}\mathrm{n}}_i\left(t_{i,coag}^{fict}\right)={\displaystyle \frac{1}{{\max }_i\left(V{C}_i\right)}}$ (6)

Within the framework of a time-driven methodology such as the MMC, the population evolution is conceptualized as a Markov process [83]. Consequently, the probability of a single agglomeration event occurring for a specific particle belonging to the ith-group during the time interval Δt is given by the following exponential expression:

${\mathrm{P}}_{\mathrm{i}}\left(\mathsf{\Delta }\mathrm{t}\right)=1-{\mathrm{e}}^{-\mathrm{V}{\mathrm{C}}_{\mathrm{i}}\mathsf{\Delta }\mathrm{t}}\approx \mathrm{V}{\mathrm{C}}_{\mathrm{i}}\mathsf{\Delta }\mathrm{t}$ (7)

To determine whether an agglomeration event occurs for the tracked particle, the Nanbu method will be utilized [84]. Initially, a random number, R₁, is generated from a uniform distribution within the interval [0,1). If this number is less than the agglomeration probability mentioned before, an agglomeration event is deemed to have taken place, in accordance with the methodology described by Zhao [85]:

$R_1\le V{C}_i\times \mathsf{\Delta }t$ (8)

The subsequent stage in the execution of an individual agglomeration event involves the determination of the agglomeration partner. Given that within the MMC framework, all agglomeration events are defined as binary collisions, only a single partner originating from a jth-group requires selection. For this procedure, the Nanbu method will be utilized again. Initially, the probability of particle i undergoing agglomeration with any other fictitious particle j is calculated as follows [85]:

$P_{ij}={\left(kwt\right)}_j\times {\mathsf{\beta }}_{ij}\times \mathsf{\Delta }t$ (9)

The relationship between the agglomeration rate and the agglomeration probability is defined as follows:

$C_i\times \mathsf{\Delta }t=\left[{\sum }_{j=1}^{N_f}\left({\mathsf{\beta }}_{ij}\times {\left(kwt\right)}_j\right)\right]\times \mathsf{\Delta }t={\sum }_{j=1}^{N_f}{P}_{ij}$ (10)

The previously determined random number R₁ can be used for the selection of the agglomeration partner. Provided that the subsequent relation is fulfilled, a particle from a jth-group will be selected as the agglomeration partner of the tracked particle of the ith-group:

${\sum }_{k=1}^{j-1}{P}_{ik}\le R_1\le {\sum }_{k=1}^j{P}_{ik},j\in \left[1,N_f\right]$ (11)

After the determination and execution of an agglomeration event, wherein two particles collide and agglomerate into a single larger one, a notable consequence is the progressive reduction of the total number of real particles: smaller particles diminish as larger ones emerge. In the context of the MMC methodology, which involves tracking fictitious particles, the numerical decrease of real particles is reflected through the appropriate adjustment of the weight of the tracked fictitious particles. As a result, the number of fictitious particles remains constant throughout the solution of the PBE, regardless of the quantity of agglomeration events considered or the temporal extent of the PBE's evolution. The operational principle employed in the present work modifies the size and number weight of exclusively one of the coagulating particles, as detailed in [86]:

${\left({\left(kwt\right)}_i\right)}_{new}={\left({\left(kwt\right)}_i\right)}_{old}\times {\displaystyle \frac{{\upsilon }_i}{{\upsilon }_i+{\upsilon }_j}}$ (12)

${\left({\upsilon }_i\right)}_{new}={\upsilon }_i+{\upsilon }_j$ (13)

${\left({\left(kwt\right)}_j\right)}_{new}={\left({\left(kwt\right)}_j\right)}_{old}$ (14)

${\left({\upsilon }_j\right)}_{new}={\upsilon }_j$ (15)

Analogously to the treatment of agglomeration events, the determination of whether a deposition event occurs for a fictitious particle i is treated as a Markov process. Consequently, the probability of its deposition within a given time interval Δt is defined by means of the deposition rate R(υ) [81]:

${{\mathrm{P}}_{\mathrm{i}}\left(\mathsf{\Delta }\mathrm{t}\right)}_{deposition}=1-{\mathrm{e}}^{-\mathrm{R}\left(\upsilon \right)\mathsf{\Delta }\mathrm{t}}\approx \mathrm{R}\left(\upsilon \right)\times \mathsf{\Delta }\mathrm{t}$ (16)

For a deposition event to be actualized, a new random number, denoted as R₂​, is drawn from a uniform distribution within the interval [0,1). The deposition will then transpire provided that the following criterion is satisfied:

$R_2\le \mathrm{R}\left(\upsilon \right)\times \mathsf{\Delta }t$ (17)

Upon the determination of a deposition event, the tracked fictitious particle must be removed from the distribution's particle population. However, this action would inevitably lead to a continuous reduction in particle number, thereby compromising the statistical treatment of the problem, a predicament identical to that encountered during the agglomeration process. To mitigate this issue, a recalibration of the weighting factor is necessary [81]. In the present work, the weight of the tracked particle is set to zero, and the simulation proceeds with this particle included as a "ghost" entity.

3. Agglomeration Kernels

The agglomeration and deposition kernels contain all the information required for modeling a process initiated by a specific physical phenomenon. Particularly for particle agglomeration, these kernels represent the probability of an agglomeration event occurring between two particles. Any mechanism that induces a relative velocity between aerosol particles has an associated agglomeration rate, and thus an agglomeration kernel [87]. The primary mechanisms driving an agglomeration process are Brownian diffusion [74], turbulence-induced motion [88,89], and gravitational settling [90], although any additional process that induces relative motion between particles, such as acoustic forcing, can also promote agglomeration. The kernel employed for solving the PBE is the total agglomeration kernel, which is composed of individual kernels representing the effect of each of the above processes. Each kernel is a function of the sizes of the two particles involved in the collision procedure. In this study, the effects of turbulent and gravitational agglomeration will not be considered since the particles examined are relatively small (micron-sized and sub-micron) and the effects of these mechanisms can be neglected [91].

3.1. Brownian Agglomeration Kernel

The movement of small particles within a quiescent fluid is characterized by the relationship of their size to the mean free path of medium molecules, as described by the Knudsen number, Kn [92]:

$Kn={\displaystyle \frac{meanfreepath}{particleradius}}$ (18)

Due to their small size, collisions between these particles and medium molecules result in their random motion (Brownian motion), leading to collisions (thus agglomeration) among the particles themselves. Depending on the value of the Knudsen number (Kn), three particle size regimes are defined: the continuum regime (Kn << 1), the transition regime (1 < Kn < 10), and the free-molecule regime (Kn > 10), and for each regime, a different expression for the collision kernel can be formulated. For the particles examined in the present work, the expression that is valid for the particles in the continuum regime will be used [93]:

${\mathsf{\beta }}_{\mathrm{c}\mathrm{o}}\left({\mathsf{\upsilon }}_1,{\mathsf{\upsilon }}_2\right)={\mathrm{K}}_{\mathrm{c}\mathrm{o}}\left({\mathsf{\upsilon }}_1^{{\displaystyle \frac{1}{3}}}+{\mathsf{\upsilon }}_2^{{\displaystyle \frac{1}{3}}}\right)\left({\displaystyle \frac{{\mathrm{C}}_{\mathrm{c}}\left({\mathsf{\upsilon }}_1\right)}{{\mathsf{\upsilon }}_1^{\frac{1}{3}}}}+{\displaystyle \frac{{\mathrm{C}}_{\mathrm{c}}\left({\mathsf{\upsilon }}_2\right)}{{\mathsf{\upsilon }}_2^{\frac{1}{3}}}}\right)$ (19)

${\mathrm{K}}_{\mathrm{c}\mathrm{o}}=2{\displaystyle \frac{{\mathrm{k}}_{\mathrm{B}}\mathrm{T}}{3\mathsf{\mu }}}$ (20)

where k_B is the Boltzmann constant, T is the absolute temperature, μ is the dynamic viscosity of the medium, and C_c(υ1), C_c(υ2) are the Cunningham slip correction factors.

3.2. Acoustic Agglomeration Kernel

The orthokinetic and hydrodynamic mechanisms have been identified as the primary contributors to the acoustic agglomeration process. To accurately simulate the impact of each mechanism on the solution of the PBE, a distinct collision kernel must be derived for each one, and the total agglomeration kernel can be expressed as the linear superposition of these kernels:

${\mathsf{\beta }}_{total}=b_{OI}{\mathsf{\beta }}^{OR}+b_{HI}{\mathsf{\beta }}^{HY}$ (21)

where b_OI and b_HI are the weighting factors of each mechanism’s contribution to the overall kernel. The multiplication factors cannot be theoretically determined. Therefore, a combination satisfying the condition b_OI+b_HI = 1 can be chosen, or alternatively, values of unity (b_OI = b_HI = 1) may be adopted, which is the approach employed in the present work.

3.2.1. Orthokinetic Agglomeration Kernel

Regarding the orthokinetic mechanism, it is crucial to analyze the oscillatory motion of particles induced by the acoustic wave. Specifically, larger particles tend to remain nearly stationary, whereas smaller particles are almost entirely entrained by the medium, oscillating back and forth with the medium [94]. The concept of the acoustic agglomeration volume was introduced by Mednikov [42] to describe the volume swept by a given particle due to its oscillatory motion (relative to other particles) during a complete acoustic cycle. Initially, large particles (collectors) are located at the center of the agglomeration volume, while small particles are distributed uniformly throughout. During a quarter of the acoustic wave period, larger particles move to the left side of the agglomeration volume relative to smaller particles, consequently capturing them. At t = 0.5T, the large particles return to the center of the agglomeration volume. During the subsequent quarter-cycle, the larger particles migrate toward the opposite side of the agglomeration volume, collecting smaller particles. At t = T, the oscillating medium conveys the larger particles back to their initial position. Mednikov’s main assumption was that all particles within the agglomeration volume are captured and that this volume is entirely replenished by other particles within a single acoustic cycle. Considering only two types of particles (1 is the oscillating collector and 2 are the particles that are inside the agglomeration volume), then the agglomeration kernel for this mechanism can be expressed as follows (Mednikov formulation):

${\beta }_{Mednikov}^{OR}={{\displaystyle \frac{1}{2}}\left(d_1+d_2\right)}^2{U}_0{n}_{12}$ (22)

where d₁ and d₂ represent the particle diameters, U₀ denotes the velocity amplitude of the sound wave, and n₁₂ is the relative entrainment factor between the two particles:

$n_{12}=H_{12}=\left|H_1-H_2\right|$ (23)

Figure 3. Schematic illustration of the acoustic agglomeration volume.

Figure 3. Schematic illustration of the acoustic agglomeration volume.

In the above equation, the H-function is the ratio of the complex velocity amplitude of the particle U_p to that of the fluid medium U₀ [95]:

$H=3{\mathsf{\lambda }}_p{\displaystyle \frac{3\left(1+y\right)-i\left(2y+3\right)y}{\left(1+y\right)9{\mathsf{\lambda }}_p-i\left[\left(2+{\mathsf{\lambda }}_p\right)2y^2+y9{\mathsf{\lambda }}_p\right]}}$ (24)

$y={\displaystyle \frac{d}{2{\delta }_v}}=\sqrt{{\displaystyle \frac{\omega {\rho }_f{d}^2}{8\mu }}}$ (25)

where λ_p is the ratio of the fluid density, ρ_f, to the particle density, ρ_p, ω is the angular frequency of the acoustic wave, and δ_v is the penetration depth of the viscous wave [96]:

${\delta }_v=\sqrt{{\displaystyle \frac{2\mu }{{\rho }_f\omega }}}$ (26)

The magnitude and phase of the H-function are respectively the entrainment factor and phase delay of the particle [97]:

${}^n=\left|H\right|=3{\mathsf{\lambda }}_p\sqrt{{\displaystyle \frac{4y^4+12y^3+18y^2+18y+9}{4{\left(2+{\mathsf{\lambda }}_p\right)}^2{y}^4+36{\mathsf{\lambda }}_p{y}^3\left(2+{\mathsf{\lambda }}_p\right)+81{\mathsf{\lambda }}_p^2\left(2y^2+2y+1\right)}}}$ (27)

$\mathsf{\varphi }={\tan }^{-1}\left[{\displaystyle \frac{\mathrm{Im}\left\{H\right\}}{\mathrm{Re}\left\{H\right\}}}\right]$ (28)

$\tan \mathsf{\varphi }={\displaystyle \frac{12\left(1+y\right)y^2\left(1-{\mathsf{\lambda }}_p\right)}{4y^4\left(2+{\mathsf{\lambda }}_p\right)+12y^3\left(1+2{\mathsf{\lambda }}_p\right)+27{\mathsf{\lambda }}_p\left(2y^2+2y+1\right)}}$ (29)

The expression for the relative entrainment factor n₁₂ can now be written as:

$n_{12}=\sqrt{n_1^2+n_2^2-2n_1{n}_2\mathrm{c}\mathrm{o}\mathrm{s}\left({\phi }_1-{\phi }_2\right)}$ (30)

A major limitation of the Mednikov model is that the oscillatory motion of particles and the mutual scattering interactions between them are not considered in the derivation of the agglomeration volume [98]. Later, Song included them in his study and modified the Mednikov expression by introducing a new formulation for the relative entrainment factor of two particles [65]:

$n_{s12}=\left|H_1-H_2+G_{12}\right|$ (31)

where G₁₂ is the so-called scattering interaction function and a detailed description can be found in Song’s dissertation.

3.2.2. Hydrodynamic Agglomeration Kernel

Apart from the linear particle-gas interactions inherent in the orthokinetic mechanism, nonlinear inter-particle interactions also occur within an oscillating gas [99]. These interactions can be attributed to either nonlinear interactions between particle scattering waves and the incident field, termed the MRPE [100], or to viscous asymmetries in the flow field surrounding the particles, leading to the AWE [101]. To investigate the influence of these interactions on particle motion, it is necessary to calculate the forces exerted on particles in the Stokes flow regime (giving rise to Bjerknes or König forces) and in the Oseen flow regime. Bjerknes forces are derived from Bernoulli's hydrodynamic principle and are applicable for arbitrary particle orientations. Conversely, Oseen hydrodynamic forces arise from the viscous asymmetry of the flow field around oscillating particles. Both force types act concurrently on particles exposed to an acoustic field; however, distinct agglomeration kernels have been formulated for each to facilitate further investigation into their impact on the acoustic agglomeration process.

3.2.2.1. Mutual Radiation Pressure Effect Kernel

The proximity of two moving particles within an acoustic field induces additional scattered waves that interfere with their trajectories. The nonlinear interaction between the incident wave and the waves scattered by each particle generates a resultant wave, disturbing nearby particles and influencing their motion. The non-uniform spatial distribution of scattering velocity leads to a net pressure force between the particles [59,102]. These forces, stemming from viscous interactions between particles and the surrounding fluid, are commonly referred to as secondary acoustic forces [103]. Such forces can result in either attraction or repulsion between proximate particles, depending on the orientation of the lines connecting the particles’ centers and the acoustic vector [91]. These forces induce a relative velocity that may culminate in collision and subsequent agglomeration [59]. In an ideal fluid, the influence of secondary acoustic forces on particle motion is considerably less significant than that of the particle's penetration into the acoustic field, leading Mednikov to disregard them [104]. Nevertheless, when particle size diminishes (consequently intensifying the effect of viscosity) or when the acoustic field is intense (at a minimum, approximately 160 dB), it becomes necessary to consider these forces and incorporate their action into the acoustic agglomeration process by formulating the corresponding agglomeration kernel.

For this hydrodynamic mechanism, the agglomeration frequency function can be described by the classic expression of König, whose calculation of the forces between two particles was based on the Bernoulli effect. König's expression for the collision kernel takes the form:

${\mathsf{\beta }}_{\mathrm{Konig}}^{\mathrm{HY}}={\displaystyle \frac{\sqrt{3}\mathsf{\pi }{\mathsf{\rho }}_f{U}_0^2}{36\mathsf{\mu }}}{\displaystyle \frac{d_1^2{d}_2^2}{d_1+d_2}}$ (32)

Danilov and Mironov [105] were the first to incorporate the viscous effects on such Bernoulli interactions and formulated an overall approach for the Bjerknes forces. Later, Song extended their work by also incorporating the effect of particle entrainment and the separation distance between particles, and formulated the following expression for the hydrodynamic mechanism [65]:

${\mathsf{\beta }}_{\mathrm{Song}}^{\mathrm{Hy}}={\displaystyle \frac{\sqrt{3}{\mathsf{\rho }}_0{U}_0^2}{144\mathsf{\pi }\mathsf{\mu }}}{\displaystyle \frac{d_1^2{d}_2^2}{d_1+d_2}}{g}_{12}=\left({\displaystyle \frac{g_{12}}{4{\mathsf{\pi }}^2}}\right){\mathsf{\beta }}_{\mathrm{Konig}}^{\mathrm{Hy}}$ (33)

The parameter g₁₂, a complex algebraic expression known as the hydrodynamic interaction function, accounts for the effects of entrainment and scattering waves; a detailed derivation can be found in Song's dissertation [65].

3.2.2.2. Acoustic Wake Effect Kernel

An asymmetry in the flow field surrounding a particle undergoing oscillatory motion in an acoustic field generates a wake characterized by an area of reduced pressure downstream of the particle. This phenomenon (acoustic wake effect) manifests under specific flow conditions, notably the Oseen flow regime, i.e., a moderate Reynolds number [106]. The presence of a proximate oscillating particle within the wake of a leading particle results in a reduction of drag, causing it to move towards the leading particle. These interacting particles exhibit mutual perturbation of their trajectories due to their periodically varying wakes [107]. Mathematical formulations for their convergent velocity in Oseen flow fields were initially developed for particles of identical size by Pshenai-Severin [101] and subsequently extended to particles of differing sizes by Dianov [108].

Though the acoustic wake effect is a purely hydrodynamic mechanism, its associated collision kernel can be expressed by considering it as a refill mechanism for the acoustic agglomeration volume, which was previously used to formulate the orthokinetic agglomeration kernel. In Mednikov's classical theory, after one acoustic cycle, all particles are swept by the core oscillating particle, and the emptied agglomeration volume is assumed to be automatically refilled by smaller particles with an arbitrarily set unity refill factor [109]. However, if physical mechanisms such as gravity, acoustic wave scattering, mutual radiation pressure, acoustic streaming, or the acoustic wake effect are considered as refill mechanisms, the volume will only be partially refilled. If the acoustic wake is considered the primary refill mechanism, the agglomeration kernel can be calculated as follows [61]:

${\beta }_{\mathrm{Hoffmann}}^{\mathrm{AWE}}={\beta }_{Mednikov}^{OR}\cdot R_{12}$ (34)

The refill factor R₁₂ has been calculated as an extension to the length of the acoustic agglomeration volume formulated by Mednikov, L₂₁:

$R_{12}={\displaystyle \frac{2L_{AWE}}{L_{21}}}=\sqrt{1+{\displaystyle \frac{6U_0\left(d_1{l}_1+d_2{l}_2\right)}{\mathsf{\pi }{\left(L_{21}\right)}^2}}{\displaystyle \frac{T_{ac}}{2}}}-1$ (35)

where l₁ and l₂ are the particle slip coefficients in the Oseen regime [108] and T_ac is the acoustic period. The length L₂₁ can be calculated from the relative displacement between the particles [109]:

$L_{21}=\left|x_1-x_2\right|={\displaystyle \frac{n_{12}{U}_0}{\mathsf{\omega }}}$ (36)

In this equation, x₁ and x₂ are the particle displacement amplitudes and L_AWE is the acoustic wake volume length that can be calculated as follows [61]:

$L_{AWE}=\sqrt{{\left({\displaystyle \frac{1}{2}}{L}_{21}\right)}^2+{\displaystyle \frac{3U_0\left(d_1{l}_1+d_2{l}_2\right)}{2\pi }}{\displaystyle \frac{T_{ac}}{2}}}-{\displaystyle \frac{1}{2}}{L}_{21}$ (37)

4. Results and Discussion

To investigate the influence of acoustic agglomeration on particle size distribution, a representative cell containing up to 4 × 10¹² real particles/m³ is isolated from an aerosol system and used in the simulation. Unless specified otherwise, the simulation is conducted with N₀ = 9,000 fictitious particles, and the size distribution is discretized into 500 logarithmically spaced bins within the range of 0.1 μm–500 μm. For the purpose of simulation validation, the MMC model will first be compared against experimental data (accounting for both Brownian motion and acoustic agglomeration) and subsequently against numerical results obtained by other researchers, based on the same experimental dataset.

The experimental data employed pertain to the temporal evolution of an initial aerosol mixture subjected to an acoustic field. These data were originally documented in the work of J. Magill [110] and further elaborated upon by P. Caperan and J. Somers [111,112], who built upon Magill's original research. These experiments investigated the time-dependent evolution of a distribution of liquid glycol particles under the influence of an acoustic field with varying intensity and frequency. All studies mentioned above utilized a standard experimental configuration. This setup comprised an agglomeration chamber for introducing the aerosol mixture, a sound source connected to a piezoelectric transducer for generating the acoustic field, and the requisite instrumentation for monitoring the particle distribution over time. The sound source is capable of operating at different frequencies and is supplied with electrical power, which can be correlated with various acoustic velocity amplitudes (as detailed in Table 1). These amplitudes were calculated using charts, provided by Caperan and Somers, that convert electrical power to acoustic pressure. Magill, in his experimental study of agglomeration, employed an acoustic chamber in the form of a 0.6 m³ acrylic glass/PMMA container. In contrast, Caperan and Somers provide a more detailed description of a slightly larger cylindrical container, measuring 2 m in length and 0.7 m in diameter. A key distinction between the two setups is that Caperan and Somers incorporated a layer of absorbent material at one end of their chamber to simulate the effect of a traveling acoustic wave on the particle mixture. Both initial mixtures are assumed to follow a log-normal distribution, with their characteristics presented in Table 1.

In addition to the experimental data, results from the work of O. Ezekoye and Y. Wibowo [59] and C. Sheng and X. Shen [113] are employed. Both studies numerically investigated the same instances of acoustic agglomeration. The former utilized a modified numerical model of Gelbard and Seinfeld’s sectional code MAEROS [114], whereas the latter solved the problem stochastically using the DSMC methodology. These findings will also serve as a foundation for validating the proper functioning of the MMC methodology in modeling acoustic agglomeration, particularly concerning the form of the agglomeration kernel, as such information cannot be gleaned from experimental results.

The liquid nature of the aerosol particles (droplets) implies that their agglomeration results in the formation of larger spherical particles of uniform density. This characteristic substantially simplifies the simulation, as all kernel expressions for spherical particles remain applicable. Furthermore, a simplifying assumption in the numerical calculations is that every collision event results in agglomeration, thus setting the collision efficiency to unity [115]. Additionally, the computational code inherently assumes spatially homogeneous conditions for both particle number concentration and the acoustic field. This necessitates a traveling wave acoustic field, a condition adopted in other simulation studies [59], and a uniform particle distribution throughout the entire field (a homogeneous mixture, facilitated by a fan within the acoustic agglomeration chamber). Lastly, as noted by Caperan and Somers, other phenomena stemming from the sound's influence on the fluid, such as acoustic streaming or acoustic turbulence, as well as the effect of the primary acoustic force, are of negligible intensity and are therefore excluded from the calculations [112].

Figure 4 compares the strengths of various agglomeration mechanisms examined in the present work. The graph pertains to a scenario where one particle has a diameter of 1 μm, the sound intensity is 100 W, and the acoustic frequencies are 10 kHz and 21 kHz. The size of the second particle ranges from 0.1 μm to 50 μm. The kernel's intensity directly correlates with its contribution to the overall agglomeration process. This figure clearly shows that acoustic mechanisms (both orthokinetic and hydrodynamic) contribute significantly more than the Brownian mechanism, with their intensity nearly two orders of magnitude higher. As anticipated, the Mednikov approach exhibits the highest intensity, which is attributable to its assumptions of an ideal fluid and complete volume refilling after a single acoustic cycle. This overestimation of the agglomeration rate will subsequently be demonstrated through comparison with experimental data. Conversely, Song’s mechanism, which incorporates the effect of scattering waves that reduce the number of particles within the agglomeration volume, is considerably weaker. Further observation of these curves reveals that the orthokinetic mechanism, in both its formulations, asymptotically approaches zero when the two particles are of identical size. This outcome aligns with its algebraic expression, as identically sized particles oscillate in the same manner within the sound field, precluding the development of relative velocity.

On the other hand, the hydrodynamic mechanism exhibits nearly constant growth with increasing particle size, and for particles of identical size, the kernel does not reach zero. This accounts for the experimentally observed agglomeration of similarly sized particles. Nevertheless, employing Song’s expression for this mechanism reveals a behavior similar to the orthokinetic mechanism, albeit less pronounced. This resemblance is anticipated, given Song's inclusion of particle entrainment in his formulation. Furthermore, with the incorporation of viscosity and the influence of scattering waves, the hydrodynamic mechanism can become comparable to the orthokinetic mechanism, potentially exerting a more significant impact on the agglomeration process.

Additionally, the expression for the kernel used to describe the acoustic wake effect is plotted in Figure 4. Since this kernel was derived as a correction to the Mednikov expression, a similar order of magnitude is anticipated. However, for equally sized particles, a reversed trend relative to the orthokinetic mechanism is observed, as the number of collisions is not diminished. This is attributed to the fact that the convergent velocity between the two particles, when calculated using the acoustic wake effect theory, is not zero for particles of the same size. From this figure, it is also evident that this mechanism is stronger than the MRPE, which is consistent with experimental observations of particle trajectories of similar size under acoustic conditions in the Oseen flow regime [116,117].

Figure 5. Contour plots illustrating the intensity of the various acoustic agglomeration mechanisms for different pairs of particle diameters. The collision frequency is expressed in units of m³/s, and acoustic velocity amplitude U₀ = 0.98 m/s with a sound frequency of (a) 10 kHz; (b) 21 kHz.

Figure 5. Contour plots illustrating the intensity of the various acoustic agglomeration mechanisms for different pairs of particle diameters. The collision frequency is expressed in units of m³/s, and acoustic velocity amplitude U₀ = 0.98 m/s with a sound frequency of (a) 10 kHz; (b) 21 kHz.

Initially, the outcomes of the MMC model, considering only Brownian motion as the agglomeration mechanism, will be presented. To monitor the temporal evolution of the agglomeration process, the non-dimensional particle number concentration, expressed as N/N₀ (where N₀ represents the initial total particle number concentration), will be depicted as a function of time. Experimental data from Caperan [112] and numerical simulation results from Sheng's DSMC method will be utilized [113]. As shown in Figure 6, both numerical simulations yield nearly identical results but underestimate the particle number concentration reduction compared to the experimental results. This discrepancy may be attributable to either the deposition of particles on the chamber walls or the gravitational sedimentation of particles. Inclusion of deposition effects in the simulation demonstrates an increased particle removal rate, aligning more closely with experimental observations. However, to ensure direct comparability with other simulation efforts, subsequent results of the acoustic agglomeration processes will exclude the influence of the deposition mechanism.

4.1. Influence of Agglomeration Kernel Expression and Residence Time

Figure 7 and Figure 8 display the temporal evolution of the non-dimensional particle number concentration, incorporating the effect of acoustic agglomeration in the simulation. A variety of individual kernels were examined, considering two distinct acoustic field frequencies at one fixed sound intensity (100 W). The combination of two kernels was performed using a simple additive method with unity multiplication factors, as previously described. In all numerical attempts, the influence of the Brownian kernel was included, but the deposition mechanism was neglected. The results from both the DSMC and MMC simulations are presented alongside the experimental data from Caperan et al. [112]. Regarding the validity of the MMC simulation, both numerical simulations predict nearly identical effects of the acoustic agglomeration process for each combination of mechanisms employed, a result that was anticipated given the nature of the two methodologies (both are statistical). A discrepancy is observed in the results of these two methodologies when the Song hydrodynamic mechanism is utilized, but this deviation is minor (4%–5% of the total range) and may be attributed to uncertainties in the formulation of the hydrodynamic interaction function employed in the DSMC. These figures illustrate that the total particle number concentration decreases over time; this decrease is significantly influenced by the acoustic field rather than by Brownian motion alone. The total number concentration is observed to decrease monotonically as the residence time increases. Nevertheless, the effect of agglomeration is more pronounced during the initial stages of the acoustic treatment compared to later periods.

The choice of kernel expressions significantly influences the simulation outcomes. For the orthokinetic mechanism, the Mednikov formulation consistently overestimates the agglomeration rate. This is because it assumes the agglomeration volume is completely refilled with particles after one acoustic cycle. In contrast, Song's approach refines Mednikov's formulation by incorporating both scattering effects and the oscillatory motion of the particle. These factors reduce the number of particles that refill the agglomeration volume, leading to fewer inter-particle collisions and, consequently, a smaller decrease in the particle number concentration. Regarding the hydrodynamic mechanism, the König formulation underestimates the experimentally observed reduction in particle number concentration, because the ideal fluid assumption neglects viscous effects, which, according to Song, could lead to an increase in the strength of the agglomeration kernel by up to eight orders of magnitude [65]. When these effects are incorporated (via Song's expression), the reduction in particle number concentration is overestimated, particularly during the initial stages of the simulation for both methodologies. However, within the context of the MMC simulation, a better agreement with the experimental results is observed during the later stages of the agglomeration process, especially for the low-frequency case.

Should the kernel expression based on the AWE be employed for the simulation, the reduction in particle number concentration is comparable to, though smaller than, the one predicted using the Mednikov approach. This is expected, as this formulation is derived from the Mednikov approach but incorporates a refill mechanism for the agglomeration volume that does not take unity values. For both cases investigated (10 kHz and 21 kHz), the results exhibit closer agreement with the experimental data. Furthermore, for the higher frequency case, the deviation from the Mednikov simulated results is less pronounced than for the lower frequency case.

Figure 9 compares the particle size distribution after acoustic agglomeration as predicted by the present model and the numerical simulation utilizing the DSMC. For this comparison, an acoustic treatment of 400 s was selected, and two distinct scenarios were analyzed: first, employing only Song's orthokinetic mechanism, and subsequently, incorporating the hydrodynamic mechanism via the König expression. The results confirm the formation of larger particles at the expense of smaller ones. Furthermore, the acoustic treatment shifts the particle distribution from a log-normal to a bimodal form, an observation consistent with acoustic agglomeration experiments conducted by other researchers [57]. When the hydrodynamic mechanism is included in the simulation, the distribution curve exhibits a flatter profile toward the higher end, leading to the creation of larger particles within the same acoustic treatment duration compared to the exclusive action of the orthokinetic mechanism. This outcome is attributed to the hydrodynamic mechanism promoting collisions between larger, similarly sized particles, which subsequently act as efficient capturers for smaller ones.

To further examine the temporal evolution of the particle distribution characteristics during the numerical simulation, the mean mass and number diameters of the particle size distribution are presented as a function of time in Figure 10. For both kernel configurations, the mass mean diameter consistently increases, indicating the formation of larger particles. When both kernels are utilized, the rate of increase is steeper. As previously explained, this is attributed to the hydrodynamic mechanism, which favors the agglomeration of particles of similar, larger sizes, while the orthokinetic mechanism promotes the agglomeration of particles with large size differences. This trend is also evident in the time evolution of the number mean diameter presented in Figure 10: since the newly created larger particles are fewer in number, their formation leads to a reduction in the number mean diameter from a certain time onwards.

4.2. Influence of Sound Field Intensity

Figure 11 illustrates the numerically predicted particle number concentration following acoustic agglomeration at various acoustic intensities. Greater acoustic intensity corresponds to an increased gas velocity, and consequently, more intense relative motion between the particles. As anticipated, the sound intensity directly influences the acoustic agglomeration results: a stronger sound field leads to a larger increase in acoustic agglomeration efficiency. This relationship is also evident from the relevant kernel expressions: the orthokinetic mechanism is directly proportional to the acoustic wave amplitude U₀, while the hydrodynamic mechanism exhibits a quadratic dependence. This figure employs the combined Song’s and König’s formulation for the orthokinetic and hydrodynamic mechanisms.

In Figure 12, a comparison is presented between the numerical simulation results for acoustic agglomeration and the experimental data reported by Magill and Caperan, considering Song's hydrodynamic expression in the numerical model. The convergence between the two numerical methodologies is again observed; however, a divergence from the experimental data is apparent across all examined cases, in which Song's expression for the hydrodynamic mechanism serves as the agglomeration kernel.

The divergence of the results of both numerical methodologies from the experimental data may be attributable not only to the algebraic form of the kernel but also to other factors affecting particle number concentration, such as the deposition or breakup of aggregates in high-intensity acoustic environments [118], which have not been included in the modeling. Furthermore, other phenomena that may enhance the agglomeration process but have not been modeled through the previously presented kernels are of particular interest: one of these is acoustic turbulence caused by the propagation of high-intensity sound fields. This phenomenon can also trigger particle collisions, thereby increasing the efficiency of acoustic agglomeration [119]. It has been observed experimentally that for sound field intensities exceeding 158 dB, this kind of turbulence appears, while for intensities greater than 160 dB, even shock waves can be generated in the medium when the sound is propagating [120]. Figure 13 presents the results obtained when utilizing the acoustic wake kernel as the exclusive agglomeration mechanism. For both simulated sound frequencies, the influence of the acoustic field's intensity is observed to be consistent with the previously studied mechanisms; specifically, a higher intensity correlates with a greater reduction in the particle number concentration.

4.3. Influence of Initial Particle Number Concentration and Shape of Distribution

The initial particle number concentration, N₀, constitutes a fundamental parameter influencing acoustic agglomeration, as is the case for all agglomeration processes. A higher initial particle number density serves to accelerate the rate of agglomeration. Regarding the orthokinetic mechanism, the reduction in initial inter-particle distance (resulting from increased concentration) enhances its effectiveness, as a greater number of particles are positioned within the agglomeration volume. Likewise, the hydrodynamic mechanism is intensified; since it is inversely proportional to the distance between particle centers, a reduction in this distance leads to a more pronounced interaction. Figure 14 shows the temporal evolution of particle number concentration within a specific acoustic field across six distinct cases of initial concentration. As evidenced by the slope of the individual curves, the agglomeration process intensifies as the initial particle count increases, a trend validated by both experimental data and corresponding numerical simulations. This behavior is consistent with experimental studies on acoustic agglomeration conducted by other researchers [31,110,121]. However, other experimental records indicate that beyond a certain threshold, an increase in particle number concentration may lead to a regime in which the agglomeration process is attenuated. This occurs due to the absorption of the acoustic field intensity by the particle mixture itself [122], a phenomenon that is not accounted for in the present study.

The initial shape of the particle size distribution also influences the acoustic agglomeration process. For this reason, numerical simulations were conducted for two additional different forms of the initial distribution: in the first case, the mean diameter was increased (keeping the number of particles and standard deviation fixed), while in the second case, the geometric standard deviation was varied (keeping the mean value and number of particles the same). This influence is best illustrated by examining the intersection of the particle cumulative size distribution curve and the acoustic particle entrainment coefficient curve, for two different acoustic frequencies, as shown in Figure 15. When the distribution is shifted toward larger particles, it is evident that more particles will acquire high relative velocities, which is expected to favor the agglomeration process insofar as the orthokinetic mechanism is concerned. On the other hand, a distribution with a large initial standard deviation spans a wider range of particle sizes, since the distribution curve (bell curve) is broader. Regarding the action of the hydrodynamic mechanism, this is also favored by the presence of larger particles, as the resulting hydrodynamic force is proportional to their size.

Figure 16 shows the results of the numerical simulation when the combination of Song’s orthokinetic formulation and König’s hydrodynamic formulation is used. The acoustic treatment was considered for a duration of 50 s, incorporating two distinct sound frequencies. It is evident that an increase in the mean diameter of the initial particle size distribution exerts a more substantial influence on the acoustic agglomeration than an increase in the geometric standard deviation.

4.4. Influence of Particle Density

Particle density is also a critical factor influencing acoustic agglomeration, as its value dictates the behavior within the acoustic field. Specifically, heavier particles (with a longer relaxation time) exhibit a slower response to fluctuations in the acoustic field, leading to a reduction in their entrainment coefficient. By plotting the entrainment coefficient for two distinct particle densities across two different acoustic frequencies (as illustrated in Figure 17), it is observed — again through the intersection of these curves — that heavier particles result in a larger population acquiring high relative velocities. Increased particle density enhances acoustic agglomeration, primarily through the orthokinetic mechanism. Conversely, the hydrodynamic mechanism remains unaffected by particle density, being solely dependent on particle diameter. This is because the driving force for the hydrodynamic mechanism, namely secondary acoustic forces, is exclusively governed by the asymmetry of the flow field surrounding the particle. Therefore, the decrease in particle number concentration is more pronounced for these heavier particles than for lighter ones, particularly during the later stages of the agglomeration process where larger particles are formed.

The predicted size distributions of particles of two different densities after 100 s of acoustic treatment are presented in Figure 18. It seems that for the heavier particles, the number reduction is more pronounced (in the range 0.1 μm–2 μm) for both acoustic frequencies examined.

4.5. Influence of Sound Frequency

Frequency is also a critical parameter influencing both the orthokinetic and hydrodynamic mechanisms. The preceding experimental results and numerical simulations indicate that a higher acoustic field frequency promotes agglomeration of the aerosol mixture. However, the influence of frequency on the evolution of the particle size distribution remains unclear upon initial consideration. Acoustic agglomeration investigations have spanned a broad spectrum of acoustic frequencies, ranging from the low values of 44 Hz [31] and 3 kHz [121], to the high-frequency fields of 10 kHz and 21 kHz [112,123]. While low-frequency sound sources were historically favored due to their cost-effectiveness and reliability, substantial advancements have recently been made in the development of devices capable of producing high-frequency and high-intensity acoustic fields. These innovations are primarily relevant for the industrial application of acoustic agglomeration [98].

To gain a deeper understanding of how frequency influences the action of each mechanism, one must examine its effect on the underlying driving force. It should be noted that, as demonstrated by the relevant equation, König's expression for the hydrodynamic kernel is solely dependent on the acoustic velocity amplitude and not on the sound frequency. This observation stems from the assumption that the fluid is considered ideal and that the particles are significantly smaller than the sound field wavelength. Only when accounting for viscous effects (i.e., wave scattering) and particle entrainment, as incorporated by Song into his expression for the hydrodynamic mechanism, does the frequency exert an influence on the corresponding collision kernel. Figure 19 illustrates the effect of frequency on the acoustic entrainment factor and the hydrodynamic interaction function g₁₂. From this figure, it is observed that the entrainment factor of a particle of a given diameter decreases with increasing sound frequency, indicating a reduction in the particle vibration speed. However, since the orthokinetic mechanism is directly related to the relative velocity of two particles, an optimal frequency exists, determined by the initial particle distribution, that promotes collisions. This is why the cumulative particle size distribution is also plotted in Figure 19. When the sound frequency is very low (~1 kHz), all particles are equally entrained by the acoustic waves (unity entrainment coefficient), resulting in minimal relative motion and thus few collisions. Conversely, when a higher frequency is selected (>20 kHz), a range of particle sizes exists where the entrainment factor transitions from unity (fully entrained) to a lower value. Particles of different sizes within this range will experience varying entrainment rates, thereby enhancing the probability of agglomeration.

Regarding the influence of the sound frequency on the effect of the scattering waves on the hydrodynamic mechanism, the most important property of this interaction between two particles is that these hydrodynamic forces depend on the source strength of these waves, which, in turn, is dependent on the product of (1−H₁) and (1−H₂) instead of the difference (H₁−H₂) [98,124]. Consequently, if H₁ and H₂ are smaller in magnitude, the hydrodynamic forces will be larger. This indicates that a higher sound frequency corresponds to stronger hydrodynamic interaction. Conversely, for very small particles or low frequencies, the entrainment factor approaches unity, resulting in these strengths being close to zero and the scattered waves being weak [124]. Furthermore, in the case of two particles of identical size, the term simplifies to (1−H)², which remains non-zero as long as the particles are not fully entrained. Thus, higher frequencies will also favor the development of relative motion based on the radiation pressure interaction and promote particle collisions. Additionally, higher frequencies result in a greater flow-around velocity, which enhances the distortion of the flow field and consequently increases the agglomeration kernel of hydrodynamic interaction [125].

Regarding the influence of sound frequency on the AWE (as investigated by Hoffmann [61]), the extended length of the agglomeration volume increases with higher frequencies, asymptotically approaching a constant value for any given pair of particles. Figure 20 illustrates the effect of sound frequency on the AWE kernel, alongside the corresponding Mednikov formulation. Since the AWE was conceptualized by Hoffmann as a variation of the Mednikov expression, a similar trend is observed when considering a pair of particles, where one particle has a fixed diameter of 1 μm and the other ranges from 0.1 μm to 10 μm. Fundamentally, increasing the sound frequency consistently enhances agglomeration [126].

Figure 21 illustrates the numerical simulation results regarding particle number concentration and the predicted particle size distribution after 100 s of acoustic treatment across four distinct acoustic frequencies. The utilized kernel combines Song's orthokinetic mechanism with König's hydrodynamic mechanism. Higher frequencies lead to a reduction in particle number concentration and broader distribution shapes, indicating the formation of larger particles at the expense of smaller ones.

5. Conclusions

The present study employs the Multi-Monte Carlo (MMC) method to solve the population balance equation (PBE) while accounting for particle acoustic agglomeration. The MMC technique utilizes the concept of "weighted fictitious particles" whose number is significantly smaller than the actual number of particles. The MMC method has been successfully applied to simulate various acoustic treatment scenarios using different agglomeration kernels. The agglomeration kernel incorporated models for orthokinetic and hydrodynamic mechanisms, as well as Brownian agglomeration and particle deposition. Furthermore, the contribution of the acoustic wake effect to the agglomeration process was evaluated. The results obtained from the MMC method exhibit good agreement with other numerical approaches such as Direct Simulation Monte Carlo (DSMC) and the sectional method (MAEROS), and available experimental results, thereby validating the computational precision and reliability of the MMC approach. A parametric examination investigated the influence of key parameters, such as sound frequency and intensity, along with initial particle size distribution characteristics, with the objective of numerically simulating the acoustic treatment of liquid aerosols. Sound frequency was identified as a critical parameter in controlling the particle size distribution and the number concentration of the treated aerosol. However, the optimal selection of this value requires joint consideration of the initial particle size distribution.

In summary, acoustically driven agglomeration leads to a diminished number density of fine particles and the subsequent emergence of larger counterparts. The proposed MMC methodology provides a rapid simulation tool for understanding and optimizing the acoustic agglomeration process.