Numerical Investigation of Aerosolization in the Venturi Dustiness Tester: Aerodynamics of a Particle on a Hill

Aerosolization of particles from surfaces [1] is of particular interest to the study of the dustiness of powders [2,3]. Natural examples of this phenomenon are pollen grains or dust suspended in air, and sand particles carried away by the wind [4–7]. Airborne organic (e.g., soot, pollens, molds, bacteria, virus) and inorganic (e.g., silica, asbestos) dust may adversely affect health. Successful delivery of powdered pharmaceuticals presupposes efficient aerosolization of the powder [8]. Previous work [9,10] studied aerodynamics of particles located on flat surfaces.

The National Institute for Occupational Safety and Health (NIOSH) has a comprehensive program to study the dustiness of powders [2,3]. NIOSH and other laboratories [11] extensively utilize the Venturi dustiness tester (VDT) [12], which permits dustiness measurement at high Reynolds number of small quantities of expensive or potentially dangerous powders (e.g., pharmaceuticals, nanoparticles). In the VDT configuration [13], a small quantity (mg) of powder is placed in a slender tube; a high velocity air stream is forced through the tube and over the powder hill; individual particles are lifted from the powder hill and are aerosolized; the aerosol is swept into a containment chamber where it is then sampled. Dubey et al. [13] modeled the air flows in the VDT, assuming the aerosolization process to have already taken place. This work addresses the aerosolization process itself, i.e., the lift and drag of particles from the powder hill. Our numerical approach is restricted to flow in the laminar regime. Sharma et al. [14,15] have extended this work to higher Re, where the dynamics of vortex shedding from the hill becomes important.

Many factors contribute to the detachment of particles from surfaces [16]. When a fully developed laminar boundary layer flows parallel to a substrate, the shear flow at the surface generates aerodynamic forces and moments. For the particle to detach, these aerodynamic forces must exceed the surface (e.g., capillary, adhesion) or body (e.g., gravity, van der Waals) forces which hold the particle onto the surface. The experimental measurement of these forces on isolated particles is challenging [17–24]. Numerical modeling is appropriate to study such problems.

For fine and ultrafine particles (micron and nanometer diameter), the particle Reynolds number (based on the local flow velocity and particle radius) is typically less than unity. The Stokes drag (creeping flow) for a sphere in a uniform steady incompressible flow in the absence of any wall is F_drag = 3 π μ u D, where u is the uniform freestream flow velocity, D is the sphere diameter, and μ is the fluid viscosity. F_lift = 0. O'Neill [25] found that the presence of a wall increases the drag by a factor f = 1.7009. Leighton and Acrivos [26] included the first inertial correction to creeping flow and found F_lift/F_drag = 0.287 Re. This differs from the earlier Saffman [27] result, also a lowest inertial correction to creeping flow, valid for large shear flow and no wall; the corresponding 2D analog, namely, of a cylinder in linear shear flow, was studied by Bretherton [28]. Cherukat et al. [29] discuss the crossover between the Saffman and Leighton–Acrivos regimes. The effect of a finite gap between particle and wall for 2 < Re < 250 has been studied [30–32]. Lee and Balachandar [33] studied the modification of lift and drag due to motion of the particle parallel to the boundary. The Finlay group [34–36] examined the related problem of a sphere or cylinder attached to a wall in a Blasius boundary layer, where additional lift derives from the vertical velocity near the plate as the boundary layer develops.

Multiple laminar experiments have been conducted [37–43]. The larger lift on an attached particle in a turbulent boundary layer has been studied [44,45]. The importance of aerodynamic forces on the migration and resuspension of small particles has been extensively studied [46–57].

The primary focus of the earlier studies has been an understanding of the wall-induced forces on a particle. More relevant for the VDT [2,12,13] is a particle of a finely divided solid (a powder) situated in proximity to a collection of other similar particles, namely, on a powder hill. The aerodynamics of a particle when it is located on a hill has not previously been studied.

In the VDT, powder is preloaded at the base of a “tee”-shaped tube (Fig. 1), which is open to the atmosphere (at top and right) and which is inserted (towards the left) into a cubical dustiness sampling chamber. In a short (Δt = 1.5 s) dispersion phase, the pressure in the dustiness chamber is lowered, and an air stream is induced to flow over the powder hill. This air stream (Re ∼ 2 × 10⁴) aerosolizes the powder particles from the hill, and the resulting aerosol enters the dustiness chamber (at left) as a jet. The earlier modeling study of the VDT [13] presupposed aerosolization—it considered an aerosol stream entering the dustiness chamber as a jet. In this work, we address the aerosolization process itself.

As discussed later, in Sec. 2.2, the air flow in the VDT tube is turbulent (Re_tube ∼ 2 × 10⁴). At these air flows, we expect vortices to be shed from the powder hill obstruction [58–62]; this is seen in the work of Sharma et al. [14,15]. Vortex shedding is suppressed due to the confined geometry [63–66]. Additional complications arise due to the porosity of the obstruction [67] and due to the dynamic disappearance of the obstacle (as particles are aerosolized from the hill). In the current, over-simplified, study, we restrict this investigation to low air flows (Re_sphere < 10²) so that the flows may be treated as steady-state. We seek to understand whether the usual drag and lift of the powder particle are significantly altered from their planar values due to the obstruction geometry. In particular, we seek to determine whether it is any easier to aerosolize a particle from the hill than from a flat surface.

It is important to determine if the lift and drag forces, which contribute to the particle detachment process, and the torque, which determines the particle rotation upon detachment, are modified for a nonplanar substrate. Hence, this study considers flow over a particle attached to a hill and compares it to the idealized case of the isolated particle attached to the plane wall. This is not a true many-body study but represents a first attempt to take into account the effect of a nonplanar substrate on particle aerosolization. We note that effects of a nonplanar substrate will be important for any treatment of particle detachment from rough surfaces, which, at the microscopic scale, consist of protuberances and valleys. Length limitations restrict the discussion in this paper to lift and drag; our simulation results for torque are contained in Ref. [9].

This paper is organized as follows: Section 1 introduces dustiness, the VDT, which motivates the problem of aerosolization of particles from a hill. Section 2 discusses the mathematical formulation underlying our simulations. We also discuss (Sec. 2.1,2) the peculiar geometry of a curved surface in contact with a flat surface, which necessitates displacing the particle slightly above the substrate. Grid refinement is discussed in Sec. 2.1.6. Conditions for the relevance of our simulations to the VDT are presented in Sec. 2.2. Section 3 presents the results of these simulations: 2D (cylinder on a wedge) in Sec. 3.1 and 3D (sphere on a cone) in Sec. 3.4. A major result of our study is that the aerodynamic lift on the particle is significantly enhanced at the summit of the hill. The origin of this effect is discussed in Sec. 3.3 and further in the Sec. S4 available in the Supplemental Materials on the ASME Digital Collection. Section 4 presents our conclusions and recommendations for further work on this problem.

This work studies the aerodynamic lift and drag forces on a particle positioned on a conical hill and slightly elevated above the hill surface. The particle is approximated as a sphere (in 3D) or as a cylinder (in 2D), with the cylinder orientated perpendicular to the freestream direction. All the simulations are performed in the laminar regime 1 < Re_channel < 2000, corresponding to 0.07 < Re_hill < 133 and to 0.01 < Re_cyl/sph < 13. We assume incompressible, steady flow, and we neglect the effects of additional surface and body forces. The parameters of standard temperature and pressure (STP) air are density ρ = 1.225 kg/m³ and dynamic viscosity μ = 1.78 × 10⁻⁵kg/m-sec; lengths are measured in mm.

where u_i is the fluid velocity.

Placing the particle in contact with the plane wall introduces a singular point at the contact point. This is numerically challenging for grid generation [35]. Because of this singular point, an improperly constructed grid may contain poor quality cells (high skewness, poor orthogonality, high aspect ratios). Potential solutions to this problem are to either slightly elevate the particle above the surface or slightly embed the particle into the surface. A brief discussion of this problem is provided in the Supplemental Material (Secs. S1 for 2D and S2 for 3D) on the ASME Digital Collection. In this study, we have taken the approach to elevate the test particle slightly above the wedge/conical substrate (perpendicular distance to the cylinder/sphere center x/R = 1.05).

A cylinder, oriented perpendicular to the flow, is positioned (Fig. 2) at different locations, θ, on a wedge submerged in a laminar fully developed flow. The diameter of the cylinder is D; the wedge height is L. The channel height is 15 L (so there is minimal constriction of the flow by the wedge), and the inlet/outlet boundaries are placed 15 L/20 L upstream/downstream, respectively, of the wedge center.

We have performed flow simulations and determined the lift, drag, and torque for the following cases: (a) three cylinder diameters, D = 2 mm, 1 mm, and 0.5 mm, at fixed particle location, θ, wedge half-angle, α, and height, L; (b) three wedge heights, L = 6.5 mm, 13 mm, and 19 mm, at fixed cylinder diameter, D, particle location, θ, and wedge half-angle, α; (c) three wedge half-angles, α = 30 deg, 45 deg, and 60 deg, at fixed cylinder diameter D, wedge height, L, and particle location, θ.

A sphere is positioned at different locations on a cone, submerged in a laminar, fully developed flow. The boundaries are placed as in the 2D case (Sec. 2.1.2.2) so that the upstream effects of the cone and sphere have minimal effects on the solution. The lift and drag forces and torque acting on the sphere are computed for various locations of the sphere positioned axially (θ) on and azimuthally (φ) around the cone (the latter at fixed θ).

The commercial code used in this study is ansysfluent v.16, which employs a finite volume method, where the conservation laws are used in their integral form. The incompressible continuity and momentum equations are solved by the semi-implicit method for pressure-linked equations (SIMPLE) algorithm, which iteratively solves for the coupled velocity and pressure fields [68]. The residual for the convergence criterion was set at 10⁻⁶. The second-order accurate implicit upwind and second-order implicit central difference schemes are used to discretize the convection and diffusion terms, respectively.

• No slip boundary at the particle surface, on top and bottom plates, and on the surface of the wedge/cone.

• Slip boundary on the side walls.

• Fully developed parabolic velocity profile is initialized at the channel inlet.

• Pressure boundary condition (set to atmospheric pressure) at the channel outlet.

Ansys Mesh Modular is used to generate structured grids for the flow geometries considered. Care has been taken to resolve wall normal velocity gradients in the boundary layer on the no-slip surfaces.

A detailed grid refinement analysis [69] is carried out on successively refined grids to estimate the error in the numerical results due to the discretization. The grid refinement study is constructed for the case with cylinder diameter D = 2 mm; wedge height L = 13 mm; wedge half-angle α = 60 deg; the cylinder is located at θ = 60 deg on the wedge. The area of the domain is A = 8.1737 × 10⁴ mm². The refinement ratio is 1.5. For the medium grid, there are 36 cells around the cylinder and six cells in the gap. The absolute drag and lift forces (not normalized to their upstream values) acting on the cylinder are computed for successive grids, and the percent changes are presented in Table 1. As the changes in these forces with mesh refinement are less than 1%, the coarsest grid was used for the reported parametric study.

A detailed grid refinement analysis [69] is carried out on successively refined grids to estimate the error in the numerical results due to the discretization. Sphere diameter D = 2 mm; cone height L = 13 mm; cone half-angle α = 60 deg; the sphere is located at θ = 80 deg on the cone. The volume of the domain is V = 1.5375 × 10⁷ mm³. Two refinements are studied: r = 1.27 (for the coarse-to-medium grid refinement); r = 1.58 (for the medium-to-fine grid refinement). For the medium grid, there are 36 cells around the sphere and six cells in the gap. The absolute drag and lift forces (not normalized to their upstream values) acting on the sphere are computed for these successively refined grids, and the percent changes are presented in Table 2. As the changes in these forces with mesh refinement are less than 1%, the coarsest grid was used for the reported studies.

The diameter of the VDT inlet nozzle is D_tube = 0.44 cm. During injection, the volumetric flow rate, Q = 60.0 L/min, corresponds to an average flow velocity u_av ∼ 65.8 m/s and Reynolds number Re_tube ∼ 19,900. While this is clearly turbulent flow throughout the bulk of the nozzle, a particle near the wall (i.e., within the viscous sublayer), experiences a much slower flow. Using a flow profile u ∼ u_max (1 − r/R)^1/n, with n ∼ 6 (appropriate for Re ∼ 19,900), u_max/u_av = (n + 1)(2n + 1)/2n², whence u_max ∼ 83.2 m/sec. For a d ∼ 1 μm particle resting on the tube surface, the velocity experienced at its center is u_particle ∼ u_max (d/D)^1/n ∼ 23.1 m/sec. The ratio of particle to tube Reynolds numbers is Re_particle/Re_tube = (u_particle/u_av)*(d/D) ∼ 8 × 10⁻⁵, whence Re_particle ∼ 1.6. Modeling the drag, lift, and torque with a local laminar flow is thus a reasonable first approach to this problem.

The cylinder (D = 2 mm) is positioned at various locations, θ, on the wedge (L = 13 mm, wedge half-angle α = 30 deg), as discussed in Sec. 2.1.2.2. Simulations are conducted for channel Reynolds numbers Re_channel < 1920, for which the flow remains laminar, consistent with pipe flow transition at Re_c ∼ 2300 [70].

Figure 3(a) displays the normalized drag, F_drag(θ)/F_drag(0), and Fig. 4(a) displays the normalized lift F_lift(θ)/F_lift(0) normalized by their upstream values (i.e., at θ = 0).

The normalized drag increases as the cylinder is located higher up on the wedge, attaining a maximum drag just upstream of the wedge apex (θ ∼ 90 deg). For 120 deg < θ < 150 deg, the cylinder is “shielded” by the wedge and experiences no drag. The maximum normalized drag is in the range 5–7.5.

On the other hand, the normalized lift varies dramatically. As the test cylinder is positioned at the various locations, θ, on the wedge surface, the normalized lift increases as the test cylinder moves up the upstream side of the wedge, reaching a maximum at θ ∼ 80 deg, and then decreases on the downstream side of the wedge; there is an “undershoot” (the lift is minimum and even negative for θ ∼ 100 deg) and then a “rebound” as the lift increases further downstream of the wedge. The enhancement effect is not small (in the range 8–19) and increases with increasing Re. The undershoot effect is similarly not small (in the range −5 to −16) but decreases in magnitude with increasing Re. Recall that a negative lift means that the airflow actually presses the test particle against the surface.

The maximum normalized lift is 3–4 times larger than the maximum normalized drag. Finally, while the lift experiences an undershoot (negative lift at θ ∼ 100 deg), the normalized drag exhibits no such anomaly.

We now discuss the effect of varying the wedge half-angle, α. The effects of varying the wedge height, L, and the cylinder diameter, D, are discussed in the Sec. S3 available in the Supplemental Materials on the ASME Digital Collection. For the wedge height study, the normalized drag (Fig. S5 available in the Supplemental Materials on the ASME Digital Collection) and normalized lift (Fig. S6 available in the Supplemental Materials on the ASME Digital Collection) are plotted as a function of channel Re.

Figures 3 and 4 present results for two additional wedge half-angles α = 45 deg and 60 deg, also at D = 2 mm and L = 13 mm. Both the normalized lift and drag appear to be relatively insensitive to the magnitude of the wedge half-angle, α. The increase in normalized drag, as the test cylinder is positioned at higher locations on the wedge, is preserved for these sharper wedges (Fig. 3). The general structure in the normalized lift (maximum lift at θ ∼ 80 deg, with minimum, negative lift at θ ∼ 100 deg) is preserved (Fig. 4).

The origin of the enhanced drag and lift forces, as the cylinder is positioned at different locations on the wedge, is examined in the Sec. S4 available in the Supplemental Materials on the ASME Digital Collection, for the case of Re = 1. Both forces result from cancelations in the pressure and shear contributions from different regions of the cylinder surface; these cancelations become less complete at higher elevations on the wedge.

Figure 5 shows that the variation in the normalized drag closely follows the variation in the normalized velocity parallel to the substrate. Similarly, the variation in the normalized lift closely follows the variation in the normalized velocity normal to the substrate (Fig. 6). While the normalized drag is amplified over the normalized parallel velocity (the dotted line is above the solid line in Fig. 5), the normalized lift is suppressed below the normalized perpendicular velocity (the dotted line is below the solid line in Fig. 6).

When the fully developed flow encounters the wedge base, a new boundary layer begins to develop as the flow turns up the wedge; the incompressible fluid is accelerated (flow constriction), with a concomitant upward vertical acceleration. Similarly, when the newly established fully developed flow (going up the wedge) encounters the apex of the wedge (flow expansion), the flow is decelerated, with a concomitant downward vertical acceleration.

The flat substrate is augmented by a cone, with half-angle α = 60 deg; the simulation flow remains laminar for Re_channel < 1920, consistent with pipe transition to turbulence at Re_c ∼ 2300 [70]. Figures 7 and 8 display the normalized drag and lift, where these are normalized using the lift and drag forces experienced by the sphere at θ = 0 deg.

The normalized drag profile (Fig. 7) for the flow over a sphere on the cone, as function of location θ, is very similar to the 2D results (Fig. 3(a)). As in 2D, maximal drag is observed at θ = 80 deg for all Reynolds numbers. The drag increases as the sphere is moved from the bottom to the peak on the upstream side of the cone and then decreases gradually to minimal drag in the “sheltered” location downstream (θ = 180 deg). For Re_channel = 1920, the maximum normalized drag is 4.5, which is smaller than the maximum normalized drag (7.4) in the 2D simulations.

The normalized lift profile (Fig. 8) for the sphere placed on the cone is similar to the 2D results (Fig. 4(a)). A maximum normalized lift occurs (θ = 80 deg) for the sphere on the upstream side of the cone apex; a minimum negative lift occurs for the sphere located on the downstream side of the cone (θ = 100 deg). This suggests that the sphere located near the top of the cone can be more easily lifted by the fluid, compared to the identical sphere located near the bottom of the cone. The variation in normalized lift is considerably larger in 3D than in 2D, e.g., maximum normalized lift (for Re = 1920) is 104 in 3D but only 28 in 2D.

The behavior near the summit (θ = 80 deg, 100 deg) is examined as a function of channel Re. The normalized drag is shown in Fig. 9; the normalized drag is shown in Fig. 10. In both cases, the Re variation is less pronounced than the variation with location on the hill.

We have also studied the azimuthal variation of the lift and drag forces on the sphere; we considered 6 azimuthal angles, φ, circumambulating the cone at θ = 60 deg.

Results for the normalized drag for flow over a sphere on a cone as function of azimuthal angle, φ, at θ = 60 deg, are plotted in Fig. 11. The reference drag used for normalization is the drag experienced by the sphere located at θ = 0 deg. There is only gradual variation in the normalized drag as the sphere moves azimuthally around the cone, with variation in the range 2.4–3.9

We have also studied the behavior of the normalized lift for the flow over a sphere placed at various azimuthal angles around the cone for θ = 60 deg (Fig. 12). The maximum lift occurs for the sphere on the upstream side of the cone (φ = 0 deg), and the minimum occurs for the sphere located on the downstream side of the cone (φ = 180 deg). The normalized lift varies significantly when compared to the variations in the normalized drag. The net upwards force experienced by the sphere decreases as the sphere is moved azimuthally around the cone. This suggests that a particle located on the upstream side of the hill is more easily lifted by the fluid, compared to a similar particle located at other azimuthal positions on the hill.

In order for the air flow to detach a particle from a flat surface, or from a hill, the aerodynamic lift force must exceed any adhesive forces (which considerably exceed the gravitational weight for fine and ultrafine particles). Measurement of particle adhesive forces is notoriously difficult, and the literature is extensive. For the following discussion, we consider two elegant experiments: (i) use of an atomic force microscope (AFM) tip to effect detachment [71]; (ii) detaching small particles via a Huntington bar shock pulse [72]. We look for data relevant to our simulations of silica particles.

Jones et al. [71] measured the adhesive force of three relevant systems: (i) small glass sphere (r = 20 μm) on the following flat plates: (a) hydrophilic glass: 4 μN < F_adh < 7 μN; (b) hydrophilic Si: 4 μN < F_adh < 9 μN; (ii) large glass sphere (r = 100 μm) on hydrophilic glass: 0.5 μN < F_adh < 2 μN. The ranges were the result of variations in the relative humidity. Wanka et al. [72] measured the adhesive force for silica spheres on a hydrophilic-treated polystyrene substrate: (i) d = 4 μm: F_adh = 70 nN; (ii) d = 10 μm: F_adh = 150 nN.

The aerodynamic lift force is given by F_lift= C_lift* (1/2 ρ u²) * πD². The particle velocity is given by u = Re_sphereν/D, whence F_lift= (π/2) * C_lift* ρν² Re_sphere² = 0.42 nN Re_sphere², using the density and kinematic viscosity of air.

Aerodynamic detachment from a flat surface should occur for the AFM system for the ranges: (i) small glass (20 μm): (a) hydrophilic glass: 100 < Re_sphere < 130; hydrophilic Si: 100 < Re_sphere < 150; (ii) large glass (100 μm) from hydrophilic glass: 35 < Re_sphere < 70. Similarly, aerodynamic detachment from a flat surface should occur for the Huntington shock systems: (i) 4 μm silica: Re_sphere ∼ 13; 10 μm silica: Re_sphere ∼ 19. Recall (Sec. 2.2) that a micron sized particle in the boundary layer of a flat surface (using the VDT flow conditions) experiences Re_sphere ∼ 2, which would seem to be insufficient to detach the particles. However, for a particle at the summit of a hill, the aerodynamic lift is enhanced by a factor of 10²; the required detachment Re_sphere decreases by a factor of 10, which is then of the right order of magnitude for the VDT to effect detachment.

In this work, we studied the variation in aerodynamic forces (lift and drag) on a test particle placed at various positions on a hill (2D–cylinder located on a wedge; 3D–sphere located on a cone). We have studied variations in lift and drag due to position on the hill (θ location), Reynolds number of the flow, and geometry of the hill (cone/wedge half-angle α); in 3D there is also azimuthal variation.

In the geometric models for the hill, the cylinder and sphere were minimally offset above the substrate, with the offset (measured normally from the particle center) required to obviate the otherwise encountered geometric singularity at the contact point. The incoming flow considered is steady laminar fully developed flow. The aerodynamic forces computed for each location (θ position) are normalized using the forces experienced by the particle upstream of the hill (i.e., at θ = 0 deg). In both 2D and 3D, the lift force experienced by the particle, near the peak of the hill upstream (wedge or cone), is greater than the force experienced by the particle on the planar substrate. This suggests that particles located upstream and near the peak of the obstacle tend to be lifted more easily.

A detailed parametric study was performed for different wedge angles, heights, and cylinder diameters. The results suggest that, as the cylinder diameter decreases, the normalized lift increases, whereas the normalized drag is relatively insensitive to cylinder size. The normalized lift and drag both increase with wedge height. In 3D, forces are also computed for the sphere at azimuthal locations around the cone.

An immediate extension of this work would be to modify the structure of the substrate to model a rough surface, which consists, at the microscopic scale, of hills and valleys. The dislodgement of particles from the rough surface (lift and drag) can then be modeled, taking into account the variation in lift and drag for the test particle situated upstream or downstream of the protuberance. Cheng et al. [73] have reported that a significant roughness reduces the normal pull-off force to a fraction of its smooth-surface value.

Any quantitative evaluation of particle detachment must take into account the physicochemical forces (surface, van der Waals, and electrostatic) between the particle and the surface. Extensive reviews of adhesion forces are available in Refs. [18] and [20]. Microparticles must overcome adhesion forces in order to detach from a surface. The adhesion force acting on the particle surface is distributed based on the surface irregularities present on the particle surface and its contacting surface. Particle shape and composition may be more important than the size, owing to their effect on the adhesion force [74].

This work has been limited to steady laminar flow conditions. Sharma et al. [14,15] have extended these results to the turbulent regime.

Finally, the powder hill has been modeled as a monolithic cone or wedge. A realistic powder hill has nonzero porosity and permeability and might be better modeled using Darcy flow through the hill.

We thank J. Bennett, K.H. Dunn, and D.E. Evans (all of NIOSH) for helpful discussions and for a careful reading of the paper. This work was supported in part by the NIOSH Nanotechnology Research Center. This paper reports research contained in the thesis of Palakurthi [9], which was the last thesis defense that Professor Kirti Ghia adjudicated. We recall his acute yet generous constructive criticism which enhanced the quality of that thesis. We dedicate this paper to his memory.

The findings and conclusions in this paper are those of the authors and do not necessarily represent the views of the National Institute for Occupational Safety and Health. Mention of product or company name does not constitute endorsement by the Centers for Disease Control and Prevention. None of the authors has a financial relationship with a commercial entity that has an interest in the subject of this paper.