About this Journal Submit a Manuscript Table of Contents
Modelling and Simulation in Engineering
Volume 2012 (2012), Article ID 687961, 10 pages
Research Article

Numerical Simulations of the Impact and Spreading of a Particulate Drop on a Solid Substrate

1School of Mechanical Engineering, Research Center for Aircraft Parts Technology (ReCAPT), Gyeongsang National University, Gajwa-dong 900, Jinju 660-701, Republic of Korea
2Department of Chemical and Biological Engineering and Applied Rheology Center, Korea University, Seoul 136-713, Republic of Korea

Received 27 April 2012; Revised 27 September 2012; Accepted 27 September 2012

Academic Editor: Jing-song Hong

Copyright © 2012 Hyun Jun Jeong et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.


We present two-dimensional numerical simulations of the impact and spreading of a droplet containing a number of small particles on a flat solid surface, just after hitting the solid surface, to understand particle effects on spreading dynamics of a particle-laden droplet for the application to the industrial inkjet printing process. The Navier-Stokes equation is solved by a finite-element-based computational scheme that employs the level-set method for the accurate interface description between the drop fluid and air and a fictitious domain method for suspended particles to account for full hydrodynamic interaction. Focusing on the particle effect on droplet spreading and recoil behaviors, we report that suspended particles suppress the droplet oscillation and deformation, by investigating the drop deformations for various Reynolds numbers. This suppressed oscillatory behavior of the particulate droplet has been interpreted with the enhanced energy dissipation due to the presence of particles.

1. Introduction

In this work, we present numerical simulations for the impact and spreading of a particulate droplet, containing a large number of solid particles, on a flat substrate and focus on particle effects on droplet dynamics such as spreading and recoil behaviors. Understanding such a behavior are particularly important in the industrial inkjet printing technology, which provides an alternative to the conventional lithography process for the purpose of accurate particle delivery to a solid surface. The present work is an extension of the authors’ previous study on a particulate droplet deformation/spread under the capillary force [1] to the inertial spreading problem by solving the Navier-Stokes equation during the droplet impact.

Numerous studies on droplet impact and spreading on a solid surface have been reported. Fukai and coworkers [2, 3] reported experimental and theoretical studies on the deformation of a droplet—the splat radius and height—impinging on a flat surface at a high Reynolds number () in an axisymmetric moving mesh. Šikalo et al. [4] performed experiments on drop impact and spreading and investigated the effective contact line and similar experiments by Park et al. [5] using a spherical cap improved predictions for flows with low Reynolds number, low Weber number, and a large contact angle. Pasandideh-Fard et al. [6] and Gunjal et al. [7] performed numerical simulations on the droplet impact using a volume-of-fluid (VOF) method and Khatavkar [8] and Khatavkar et al. [9] employed a diffuse-interface method to study the capillary and low inertia spreading of a micro droplet on a flat and patterned surface at low Reynolds numbers. Recently, Lee and Son [10] presented droplet impact and coalescence to understand microline patterning process with the level-set method. In most previous studies, however, both the experiments and simulations focused on drop deformation with a homogeneous fluid, whose application to the particulate fluid drop is limited due to the presence of the hydrodynamic interaction between particles and the fluid. Although interest in the particle effects as they pertain to drop deformation has increased, only a few experimental studies have explored the particulate droplet behavior. Furbank and Morris [11] observed the drop formation characteristics and Nicholas [12] studied the effect of particles on the spreading behavior, finding that a particulate drop spread less than a pure drop. Recently, the authors performed 2D direct simulations of low speed capillary spreading of a particulate droplet considering hydrodynamic interaction of particles and equilibrium contact angle [1].

In this work, we develop a direct numerical simulation technique for particulate drop spreading with the impact on a solid substrate. This work is restricted to two-dimensional simulations, owing to the complexity in the droplet impact and spread phenomena of the particulate drop. We present numerical results for the evolution of the drop shape, the particle distribution and movement, and the effects of the particles on the dynamics of droplet spreading. The term “direct” is used here because hydrodynamic and interparticle interactions are fully considered at the individual particle level. For this purpose, we employ the direct simulation scheme of Hwang et al. [13] similar to the distributed Lagrangian-multiplier (DLM) method [14]. We use the second-order Adams-Bashforth/Crank-Nicholson method (AB2/CN) to solve the Navier-Stokes equation and the level-set method to track the drop-air interface, dealing with a free-surface problem in a fixed mesh. The interfacial tension is assigned via the surface stress, known as the continuous surface stress (CSS) formulation and a virtual stress tensor is introduced at the contact point to exert force toward the equilibrium contact angle. Using the proposed numerical method, we study the impact and spreading behaviors of a particulate drop on a flat surface and discuss the effect of such particles by investigating particle movements, changes in the drop spreading and oscillation behaviors, and the kinetic energy.

2. Modeling

In the present work, we consider the impact and spreading of a 2D particulate liquid drop on a solid substrate and, as described in Figure 1, an initially circular drop with radius is in contact with a solid surface. An instantaneous artificial gravity is applied for a short time period to generate the initial velocity distribution within a drop as well as in the air domain, satisfying the incompressibility condition of the entire domain. The resultant velocity at the drop center, generated in this way, is regarded as the characteristic drop impact velocity . Therefore the drop containing a number of randomly distributed rigid circular particles undergoes impact, spreading, and recoil behaviors through the inertial force, the elastic interfacial force, and the capillary force and subsequently reaches an equilibrium state. The spreading of a droplet in 2D is a free-surface problem with a moving contact point that requires a proper representation of various interfaces, such as the air-liquid and liquid-solid interfaces, to account for the interfacial tension and wetting effects [15]. In addition, particles inside the drop experience complex flow fields during this process, and particle motions are determined by hydrodynamic interaction between particles and the surrounding drop fluid.

Figure 1: A schematic representation of a particulate droplet spreading. The initial drop contacts the solid surface and the particles are distributed inside the drop. The drop spreads with inertial force and then spreads to an equilibrium state.

The problem is described in Figure 2 along with the coordinate system. The computational domain of the size is denoted by . The four boundaries of the domain are denoted by . The region occupied by a droplet is denoted by with the interface between the drop and the air. The symbol is the outward unit normal vector on the drop interface. Particles are denoted by (i=1,…,N) and is the number of particles. We use a symbol for a collective region occupied by the particles at a certain time . For a particle , ,, ,, and are used for the locations of the particle center, the translational velocity, the angular velocity, the angular rotation, and the radius, respectively. The vector is the unit vector normal to the plane.

Figure 2: The computational domain with the droplet interface and particles inside in this study.

The set of governing equations for the region occupied by a drop at a certain time is given by Equations are equations for the momentum balance, the continuity, and the constitutive relationship of a Newtonian fluid and the symbols , , , , , , , and are the drop density, the velocity, the gravity, the stress, the pressure, the identity tensor, the drop viscosity, and the rate of deformation tensor, respectively. Accepting the Eulerian approach, we replace the air with a fictitious fluid having a viscosity of one hundredth of the droplet viscosity [16]. Then we have exactly the same set of equations in the air domain (the fictitious fluid) for the momentum balance, the continuity, and the constitutive equations as those of the droplet fluid and the air region is identified by the air viscosity and the air density . Note that the subscript “” denotes the air domain. In this sense and for convenience in future use, one can introduce quantities and such that and in the drop domain and and in the air domain. This Eulerian approach facilitates the use of single mesh for the entire computations with the presence of particles and the fluid interface.

As mentioned earlier, to apply the inertial force to a drop, we employ an instantaneous artificial gravity, and it is applied to an initially stationary drop only for a short time period to accelerate the drop and to achieve a specific initial velocity . For this purpose, we replace the gravity with the instantaneous artificial gravity in the momentum equation as follows: As for the interface condition between the drop and the air, the following stress balance should be satisfied in the absence of interfacial tension gradient: The variables and denote the interfacial tension and the local mean curvature, respectively, on the interface and is the local unit outward normal vector to the interface. To deal with the interface condition (3) in the Eulerian frame, we employ the level-set method and in this description the interface is not explicitly tracked, but is conveniently described by the level-set function . The zero-level set indicates the interface and it takes a positive value within the air region and a negative value inside the drop. In fact, the level-set function is a signed distance function such that the magnitude is determined by the shortest distance from a point to the interface: that is, This level-set function is advected passively under a given flow field and therefore its evolution in the Eulerian frame can be expressed as If the level-set function is seriously distorted as time evolves, the level-set function must be redistributed to satisfy the distance function property (4) and we employ a direct geometric redistancing scheme developed by Kim and Hwang [17] for this purpose. This level-set description allows a simple treatment of the interfacial condition (3) by introducing the smeared-out interfacial stress with the delta function based on the level-set function, which is called the continuous surface stress (CSS) formulation [18]: and it can be added to the constitutive equation as an additional stress where is the total stress and (7) defines the constitutive relation over the entire domain, including the droplet/air interface and, as we will see later, the interior of suspended particles. In addition due to our unified descriptions for drop and air domain, the momentum balance for the entire domain can be written simply as

In this problem, there are three driving forces in the droplet dynamics: the impact of droplet, elasticity associated with the interfacial tension, and the capillary force. Among them, the capillary force is the microscale action on the contact point with a solid wall, where complicated microscopic phenomena are present, leading to in an equilibrium contact angle macroscopically. In the present work, only the static equilibrium is considered following the work of Bellet [19], where an additional force at the contact point is introduced, if the current contact angle is different from the static equilibrium contact angle , to generate the motion toward static equilibrium. To combine it with the CSS formulation, we artificially enforce the interfacial stress tensor near the contact point (within around an element distance) such that using the equilibrium normal vector , the normal , and tangential , vectors to the solid wall (Figure 3).

Figure 3: Modification of the normal vector at the contact point.

Now we consider the description of freely suspended particles. Following the work by Hwang et al. [13], a circular particle is modeled as a rigid ring filled with the same fluid used in the drop fluid domain. The rigid body condition is imposed on the particle boundary only and it provides a plausible description for particles in fluid in low Reynolds number limit. This fictitious domain approach for the particle description is particularly useful in this work, given that it facilitates the use of a single Eulerian mesh for unified treatments of both the free-surface with interfacial tension and the solid/liquid interaction. The rigid-body condition on the particle boundary states that the velocity of a fluid particle attached onto the particle boundary is identical to that of the particle: Underlying assumption in this description is that the inertial and gravity forces inside the rigid ring may be safely neglected, which can be satisfied in this problem as well, in spite of the presence of the droplet impact. One way to check this applicability is to measure the particle Reynolds number and the drop Reynolds number and the particle Reynolds number , are defined as where and are the diameters of the particle and drop, respectively. In a typical inkjet printing process, the drop Reynolds number is less than 10 and, considering the drop and particle size ratio , the particle Reynolds number becomes less than , which facilitates the use of the rigid-ring description. Although not presented here, we tested another scheme, with which the rigid-body constraint is applied not only on the particle boundary but also in the particle interior, and the resultant rigid-body velocities were compared with results from the rigid-ring approximation. The two schemes were found to yield identical result.

The movement of the circular particle is given by kinematics To determine the unknown particle rigid-body motion , one needs the balance conditions for the drag force and the torque on the particle boundary. In the absence of inertia and external forces or torques, particles are force-free and torque-free: We did not use an artificial particle-particle collision scheme as in Glowinski et al. [14], as particle overlap and particle/wall collisions could be avoided through the use of a relatively small time step and sufficiently refined particle boundary discretization.

3. Numerical Methods

Jeong et al. [1] presented the weak form for capillary spreading of a particle-laden droplet using (i) the standard velocity-pressure formulation for the Stokes flow, (ii) the combined weak formulation of Glowinski et al. [14] for implicit treatment of freely suspended particles, and (iii) the level-set method with the CSS formulation [17] for the interface description between droplet and fluid. In the present work, an additional treatment is necessary to incorporate the inertia and convective terms in the Navier-Stokes equation for the description of impact and subsequent spreading behaviors. For this purpose, we employ the second-order Adams-Bashforth/Crank-Nicholson time stepping scheme to solve the Navier-Stokes equation [20]. In this scheme, the convective term is treated explicitly with the second-order Adams-Bashforth (AB2) method and the viscous term is treated implicitly with the Crank-Nicholson (CN) method. This time-stepping scheme for the Navier-Stokes equation can be written in a strong form as follows: In the above equation, we introduce an extra stress tensor that includes the viscous stress and the interfacial stress .

The numerical method can be best interpreted with the time stepping procedure due to the AB2/CN time stepping for the momentum equation as below. At every time step, we take the following four step procedures.

Step 1. For given fluid velocity , we determine the level-set function by solving the evolution equation of the level-set function: To deal with the convection term, we employ the discontinuous Galerkin method in (15), in which is the outward unit normal vector on the boundary of an element , is the part of the element boundary where (inflow), and is the level-set function value in the neighboring upwind element. Although (15) is expressed in terms of the explicit Euler method, what is actually employed is the TVD-RK3 (total variance diminishing 3rd order Runge Kutta) method by applying this explicit Euler method three times for the single time step advance [17, 21].

Step 2. Given the particle velocity , the next particle positions are updated (12).

Step 3. Given , the interfacial stress is evaluated from (6) using .

Step 4. Given ,, and , the weak form for the momentum and mass balances are solved with the constraint equation for the freely suspended particles. The weak form can be stated as follows: find such that for all admissible functions . In (16), the AB2/CN scheme is incorporated within the weak form and the rigid-ring constraint (10) has been combined as well by introducing the Lagrangian multiplier defined on the boundary of each particle. The Lagrangian multipliers can be (roughly) identified as traction forces acting on the particle by the fluid [13].

We have a few remarks on this four step procedures. At the first time step , the level-set function and particle positions are set and the problem in Step 4 is solved to obtain velocities and of particles and fluid. Then, at each time step, we repeat the four step procedures to advance to the next time step. The AB2/CN cannot be implemented at the first step as the fluid velocity at the previous time step is not known and the explicit-Euler/CN has been used instead. As the level-set function evolves over time, the distance function property of (4) is no longer satisfied and the level-set function must therefore be reinitialized to restore the distance function property, because the accuracy of the distance function is crucial for the evaluation of the interfacial stress . In this work, we reinitialize the level-set function by a direct geometrical redistancing method, following the work of Kim and Hwang [17]. We search for line segments along the interface on the subdivided element boundary and introduce a separate function, the distance function, whose function value is the shortest distance from any nodal position to the line segments along the interface. Details of the geometric construction of the distance function can be found in Kim and Hwang [17].

For the discretization of the weak form, we use a uniform regular quadrilateral element with the continuous biquadratic interpolation for the velocity , a discontinuous linear interpolation for the pressure , a discontinuous biquadratic interpolation for the level-set function , and the collocation point method the Lagrangian multipliers for the rigid-ring constraint. As the discretization method in the present study is the same as one of the author’s previous work, details of the implementation technique can be found in Hwang et al. [13], Kim and Hwang [17], Jeong et al. [1]. At every time step, we solve a large matrix equation from the discretization of the weak form (16) using a direct method based on a sparse multifrontal variant of Gaussian elimination (HSL/MA41).

4. Numerical Results

In this section, we present numerical results for the impact and spreading of a particulate droplet on a flat surface. A finite-element mesh of is shown in Figure 4 along with the initial configuration of a droplet with 50 particles inside. The computational domain is of size , the initial drop radius is and the particle radius is 140 . The drop center is located initially at (5 mm, 1.9 mm) such that the bottom part of the drop is attached onto the solid surface from the beginning, as mentioned earlier. The particle fraction, denoted by , is defined as the ratio of the particle area to the total drop area and is fixed as throughout this study. In the authors’ previous work [1], we investigated the particle-laden droplet spreading by the capillary force in the absence of inertia four different particle fractions (from zero to 0.954) and found that kinetic energy dissipation and local shear rate increase with the particle fraction, which leads to slow drop spreading even in the drop spreading processes governed by the capillary force. For this reason, we have focused in this work a relatively high particle fraction to clearly identify particle effects on drop deformation in this much more complex problem with effects of droplet impact with inertia.

Figure 4: An example finite-element mesh with the initial drop interface. Particles, denoted as collocation points along the particle boundary, are distributed randomly inside the drop.

A particle is described by uniformly distributed collocation points on its boundary and in this problem ten collocation points are used to represent a particle. Material properties are roughly taken from the data of Nicolas [12]. The density and viscosity of the drop fluid are 400 kg(cm3) and 0.01 Pa·s, respectively. As for the air region, we choose the density and the viscosity such that and . The interfacial tension coefficient is set to 2 mNm, which is ten times smaller than the experimental data of Nicolas [12], due to the numerical stability problem and the equilibrium contact angle is set to . The instantaneous artificial gravity ranges from zero to 4000, which correspond to 0 to 20.389 g with the gravitational constant . The initial impact velocity is the velocity at the droplet center, which is generated by the application of for 0.0002 [sec]. In this case, the drop Reynolds number , the Capillary number , and the Weber number are in ranges from 1.45 to 22.73, from 0.66 to 161.45 and from 0.45 to 7.1, respectively. The time step used in this study is 0.0002 [sec] for  g and 0.00005 [sec] for 20.389 .

Plotted in Figure 5 are consecutive droplet deformation patterns of a particulate droplet under the capillary spreading (Figure 5(a)) and the spreading accompanied by the impact to illustrate overall behaviors of inertial spreading in comparison with the capillary spreading ( The corresponding movie files for this numerical experiment are accompanied with the paper as supplementary data: sup_movie_Figure  5(a).avi (capillary spreading) and sup_movie_Figure  5(b).avi (impact and spreading) available online at doi: 10.1155/2012/687961). (Note that the dimensions in all figures in this work are normalized by ). The dimensionless parameters in the impact problem are the drop Reynolds number , the capillary number , and the Weber number . In the case of capillary spreading, as time passes, the drop spreads slowly and monotonically toward the equilibrium configuration. In contrast, inertial spreading looks more dynamic and indeed shows oscillation behavior: it reaches maximum spreading quickly at [sec] and then recoils back to near the original shape and the maximum recoil behavior is observed at [sec]. Subsequently, the drop spreads slowly and monotonically to the equilibrium state. To recognize the particle movement clearly, we colored three distinctive particles in black originally located near the bottom surface for both cases. As seen in Figure 5 as well as in the supplementary movies, the most evident motions of particles are generated by the impinging flow toward the surface in both cases, which is mainly driven by the spreading of the interface. In case of the impact spreading, the particle movement in the early fast spreading (before recoil) is found to be minor but it is observed apparently during the second spreading process after recoil. Particle redistribution is enhanced in the impact spreading due to high kinetic energy and active spreading in this case.

Figure 5: The consecutive drop deformation shapes with particle distributions (a) for the capillary spreading and (b) the inertial impact spreading.

Next we investigate effects of particles on the droplet oscillation behaviors. Presented in Figure 6 are droplet interfaces at the instances of maximum spreading and maximum recoil for both the homogeneous fluid droplet and particulate droplet under the inertial impact spreading. Simulation conditions are taken the same as the previous inertial spreading of a particulate drop (Figure 5(b)). The homogeneous drop is observed to reach the maximum spreading at [sec] (Figure 6(a)) and the maximum recoil at  s (Figure 6(b)), whereas the particulate drop reaches the maximum spreading slightly faster at [sec] but the maximum recoil occurs later at than the homogeneous fluid drop case. In other words, the homogeneous drop spreads more and longer during the early spreading (before recoil) but also recoils faster than the particulate drop. This implies that the presence of particles reduces the mobility of the droplet spreading and recoil dynamics. Investigating the droplet interface at the maximum spreading and at the maximum recoil, one can observe that the homogeneous fluid drop spreads more and recoils more than the particulate drop, which indicates again the reduced mobility of the particulate drop.

Figure 6: The droplet interfaces at the instances of (a) the maximum spreading and (b) the maximum recoils for a homogeneous fluid drop and the particulate fluid drops.

In order to provide physical insight of droplet spreading and oscillation, we present a high Reynolds number case (Re = 22.73, We = 161.45, and Ca = 7.1) as two additional supplementary movies: one for the pure fluid case (sup_movie_high_Re_purefluid.avi) and one with the presence of particles (solid fraction = 0.245, sup_movie_high_Re_suspension.avi).

For quantitative comparison of the droplet deformation, we present in Figure 7 the spread factor and the splat height factor for three different cases discussed in the previous examples: (i) the capillary spreading of a homogenous fluid droplet; (ii) the inertial impact spreading of a homogenous fluid droplet; (iii) the inertial impact spreading of a particulate droplet. The spread factor is the droplet width normalized by the initial drop diameter and the splat height is defined as the droplet height divided by the initial droplet diameter. From Figure 7(a) for the spread factor, the homogenous droplet with impact is found to show the largest spread in the early inertial spreading phase, whereas the capillary spreading of the homogeneous drop does not show any oscillatory dynamics. The inertial spreading of the particulate droplet shows intermediate behaviors, which is consistent with the previous observation with Figure 6. The evolution of the splat height in Figure 7(b) shows monotonic decrease in case of capillary spreading of the homogeneous fluid droplet but fast drop followed by recoil back in case of the inertial impact and spreading. In the latter case, the amount of decrease in the splat height seems to be reciprocal to the particle fraction, which indicates again that particles restrain the oscillation of the drop. We remark that the evolution of the splat height in the homogeneous fluid droplet shows qualitative agreement with the data of Khatavkar [8], though not presented here.

Figure 7: The comparison of time evolutions of (a) the spread factor and (b) the splat height for capillary and inertial spreading with the presence or absence of particles.

Suppressed oscillation due to the presence of particles may be related to energy dissipation, as discussed in Jeong et al. [1]. The presence of particles increases the shear rate and thereby the shear stress near the particles, which is expressed macroscopically by the increased bulk viscosity of the particulate fluid. The local energy dissipation can be expressed as the product of the increased bulk viscosity and square of “averaged” local shear rate (or the viscosity times square of increased shear rate due to presence of particles) and the observation of kinetic energy within a droplet during spreading may yield rough estimates for the amount of the energy dissipation. In this regard, we present the average kinetic energy in time in Figure 8, which is defined as over the drop area including the interior of particles, for several Reynolds numbers with or without particles inside. (The particle velocity can be roughly identified by the average fluid velocity surrounding the particle.) One can observe completely distinct behavior of capillary spreading and the inertial spreading. In the capillary spreading, the kinetic energy increases due to the acceleration due to the capillary force, though slow, and then as the drop shape approaches the equilibrium configuration it vanishes to zero. On the other hand, there is a fast decrease of kinetic energy during the early spreading phase and then uprise of kinetic energy can be observed due to recoil of droplet. This means that a certain amount of kinetic energy is stored as the elastic energy of the interface during the first decrease in kinetic energy. Another observation is that the initial kinetic energy increases with the Reynolds number, as it should be. More importantly, one can see kinetic energy of the particulate drop decrease faster than that of the homogeneous fluid droplet in all cases, which supports our arguments of enhanced energy dissipation due to presence of particles, and the amount of uprise due to droplet recoil is again found to decrease with the particle presence. This increased energy dissipation seems to be responsible to the suppressed oscillatory behavior of the particulate droplet.

Figure 8: The averaged kinetic energy in time during spreading for various flow conditions.

5. Conclusions

In this work, we presented a 2D numerical simulation technique for a particulate droplet spreading primarily driven by the inertial impact on a flat solid substrate and investigated particle effects on the particulate droplet dynamics such as spreading and recoil behaviors. Based on a finite-element method for the Navier-Stokes equation, the level-set method has been employed to describe the droplet/air interface dynamics and a DLM-like fictitious domain method has been introduced to describe freely suspended particles within the droplet. The combination of these techniques facilitates the use of a single Eulerian mesh in the entire computation.

We report that unlike the capillary spreading a particulate droplet reaches the maximum spreading quickly from the impact and recoils back to near original shape; and then spreads monotonically to an equilibrium state. During the course of impact and spreading, particle movements are found to be mainly driven by the latter monotonic spreading. In addition, we observed that the presence of particles reduces the mobility of the droplet spreading and recoil dynamics that invokes the reduced amount of the spread factor at the maximum spreading and slower deformation during spreading and recoil. Investigating changes in the average kinetic energy in time, we hypothesized that the increased energy dissipation due to the presence of particles would be responsible for the reduced mobility of the particulate droplet.


This work was supported by the National Research Foundation (NRF) of Korea funded by the Ministry of Education, Science, and Technology (2011-0031383 and 2010-0021614).


  1. H. J. Jeong, W. R. Hwang, C. Kim, and S. J. Kim, “Numerical simulations of capillary spreading of a particle-laden droplet on a solid surface,” Journal of Materials Processing Technology, vol. 210, no. 2, pp. 297–305, 2010. View at Publisher · View at Google Scholar · View at Scopus
  2. J. Fukai, Z. Zhao, D. Poulikakos, C. M. Megaridis, and O. Miyatake, “Modeling of the deformation of a liquid droplet impinging upon a flat surface,” Physics of Fluids A, vol. 5, no. 11, pp. 2588–2599, 1992. View at Scopus
  3. J. Fukai, Y. Shiiba, T. Yamamoto et al., “Wetting effects on the spreading of a liquid droplet colliding with a flat surface: experiment and modeling,” Physics of Fluids, vol. 7, no. 2, pp. 236–247, 1995. View at Scopus
  4. S. Šikalo, M. Marengo, C. Tropea, and E. N. Ganić, “Analysis of impact of droplets on horizontal surfaces,” Experimental Thermal and Fluid Science, vol. 25, no. 7, pp. 503–510, 2002. View at Publisher · View at Google Scholar · View at Scopus
  5. H. Park, W. W. Carr, J. Zhu, and J. F. Morris, “Single drop impaction on a solid surface,” AIChE Journal, vol. 49, no. 10, pp. 2461–2471, 2003. View at Publisher · View at Google Scholar · View at Scopus
  6. M. Pasandideh-Fard, S. Chandra, and J. Mostaghimi, “A three-dimensional model of droplet impact and solidification,” International Journal of Heat and Mass Transfer, vol. 45, no. 11, pp. 2229–2242, 2002. View at Publisher · View at Google Scholar · View at Scopus
  7. P. R. Gunjal, V. V. Ranade, and R. V. Chaudhari, “Dynamics of drop impact on solid surface: Experiments and VOF simulations,” AIChE Journal, vol. 51, no. 1, pp. 59–78, 2005. View at Publisher · View at Google Scholar · View at Scopus
  8. V. Khatavkar, Capillary and low inertia spreading of a microdroplet on a solid surface [Ph.D. thesis], Eindhoven University of Technology, 2005.
  9. V. V. Khatavkar, P. D. Anderson, and H. E. H. Meijer, “Capillary spreading of a droplet in the partially wetting regime using a diffuse-interface model,” Journal of Fluid Mechanics, vol. 572, pp. 367–387, 2007. View at Publisher · View at Google Scholar · View at Scopus
  10. W. Lee and G. Son, “Numerical study of droplet impact and coalescence in a microline patterning process,” Computers and Fluids, vol. 42, no. 1, pp. 26–36, 2011. View at Publisher · View at Google Scholar · View at Scopus
  11. R. J. Furbank and J. F. Morris, “An experimental study of particle effects on drop formation,” Physics of Fluids, vol. 16, no. 5, pp. 1777–1790, 2004. View at Publisher · View at Google Scholar · View at Scopus
  12. M. Nicolas, “Spreading of a drop of neutrally buoyant suspension,” Journal of Fluid Mechanics, vol. 545, pp. 271–280, 2005. View at Publisher · View at Google Scholar · View at Scopus
  13. W. R. Hwang, M. A. Hulsen, and H. E. H. Meijer, “Direct simulation of particle suspensions in sliding bi-periodic frames,” Journal of Computational Physics, vol. 194, no. 2, pp. 742–772, 2004. View at Publisher · View at Google Scholar · View at Scopus
  14. R. Glowinski, T.-W. Pan, T. I. Hesla, and D. D. Joseph, “A distributed Lagrange multiplier/fictitious domain method for particulate flows,” International Journal of Multiphase Flow, vol. 25, no. 5, pp. 755–794, 1999. View at Publisher · View at Google Scholar · View at Scopus
  15. J. U. Brackbill, D. B. Kothe, and C. Zemach, “A continuum method for modeling surface tension,” Journal of Computational Physics, vol. 100, no. 2, pp. 335–354, 1992. View at Publisher · View at Google Scholar · View at Scopus
  16. G. A. A. V. Haagh and F. N. Van De Vosse, “Simulation of three-dimensional polymer mould filling processes using a pseudo-concentration method,” International Journal for Numerical Methods in Fluids, vol. 28, no. 9, pp. 1355–1369, 1998. View at Publisher · View at Google Scholar · View at Scopus
  17. S. J. Kim and W. R. Hwang, “Direct numerical simulations of droplet emulsions in sliding bi-periodic frames using the level-set method,” Journal of Computational Physics, vol. 225, no. 1, pp. 615–634, 2007. View at Publisher · View at Google Scholar · View at Scopus
  18. B. Lafaurie, C. Nardone, R. Scardovelli, S. Zaleski, and G. Zanetti, “Modelling merging and fragmentation in multiphase flows with SURFER,” Journal of Computational Physics, vol. 113, no. 1, pp. 134–147, 1994. View at Publisher · View at Google Scholar · View at Scopus
  19. M. Bellet, “Implementation of surface tension with wall adhesion effects in a three-dimensional finite element model for fluid flow,” Communications in Numerical Methods in Engineering, vol. 17, no. 8, pp. 563–579, 2001. View at Publisher · View at Google Scholar · View at Scopus
  20. J. Kim and P. Moin, “Application of a fractional-step method to incompressible Navier-Stokes equations,” Journal of Computational Physics, vol. 59, no. 2, pp. 308–323, 1985. View at Scopus
  21. C. W. Shu and S. Osher, “Efficient implementation of essentially non-oscillatory shock-capturing schemes,” Journal of Computational Physics, vol. 77, no. 2, pp. 439–471, 1988. View at Scopus