Spline Based Shape Prediction and Analysis of ... - ACS Publications

Article pubs.acs.org/Langmuir

Spline Based Shape Prediction and Analysis of Uniformly Rotating Sessile and Pendant Droplets Karan Jakhar, Ashesh Chattopadhyay, Atul Thakur, and Rishi Raj* Department of Mechanical Engineering, Indian Institute of Technology Patna, Patna, Bihar, India S Supporting Information *

ABSTRACT: Prediction and analysis of the shapes of liquid− vapor interface of droplets under the influence of external forces is critical for various applications. In this regard, a geometric model that can capture the macroscopic shape of the liquid−vapor interface in tandem with the subtleties near the contact line, particularly in the regime where the droplet shape deviates significantly from the idealized spherical cap geometry, is desirable. Such deviations may occur when external forces such as gravity or centrifugal dominate over the surface tension force. Here we use vector parametrized cubic spline representation for axisymmetric fluid−fluid interfaces along with a novel thermodynamic free energy minimization based heuristic to determine the shape of liquid−vapor interface of droplets. We show that the current scheme can easily predict the shapes of sessile and pendant droplets under the action of centrifugal force over a broad range of surface contact angle values and droplet sizes encountered in practical applications. Finally, we show that the cubic spline based modeling approach makes it convenient to perform the inverse analysis as well, i.e., predict interfacial properties from the shape of a droplet under the action of various types of external forces including gravity and centrifugal. We believe that this versatile modeling approach can be extended to model droplet shapes under various other external forces including electric and acoustic. In addition, the simple shape analysis approach is also promising for the development of inexpensive interfacial analysis tools such as surface tensiometers.

■

INTRODUCTION Fundamental understanding of interfacial shape and contact angle of droplets on various surfaces and under the action of external forces is important for microfluidics, biological, manufacturing, and heat transfer applications.1−4 At small volumes, i.e., Bond number Bo = ρgL2/σ ≪ 1, surface tension forces dominate and the liquid−vapor interface of a sessile droplet placed on a chemically homogeneous and smooth surface assumes spherical cap geometry where the contact angle is dictated by Young’s equation.5 Here ρ and σ are the density and the surface tension of fluid, respectively, g is the acceleration due to gravity, and L is the characteristic dimension of the droplet (L = (3V/4π)1/3, where V is the volume of the droplet). However, in practical applications, droplets may be large (Bo ≥ 1) and the liquid−vapor interface may undergo deformation such that the shape may deviate from spherical cap geometry. Deviation from spherical cap geometry is also observed in the case of large pendant droplets (Bo < 0) or when a droplet is subjected to various types of external forces including centrifugal forces. The prediction of the shape of a droplet resting on a horizontal surface which pertains to the balance of gravity and surface tension is classical, and numerical solutions are known for more than a century.6 Further studies are continuously made to estimate the geometry of the liquid−vapor interface of © 2017 American Chemical Society

a droplet on diverse surfaces and under the effect of various types of external forces through numerical solutions of the Young−Laplace equation using the finite difference method,7 finite element method,8,9 calculus of variation,10 and spectral boundary element algorithm.11−19 For example, Dimitrakopoulos et al.19 performed a detailed analysis on the interface of droplet under the effect of gravity. They presented a comprehensive discussion on the boundary element algorithm and simulated the shapes with varying Bo and contact angles. Conversely, geometric approaches involving B-splines,20,21 elliptic contours,22 metaball,23 and CAD based algorithms24,25 have also been used. Although these methods are accurate in predicting the droplet shapes, solution formulations are rather involved. Approximate analytical solutions of Young−Laplace equation are often adopted as a convenient alternative, however, are mostly valid for small range of parameters beyond which the approximation breaks down and the shape prediction results are not accurate. For example, while some analytical solutions are available for high surface tension fluids,26,27 very little literature is available on approximate analytical solutions for low surface Received: March 10, 2017 Revised: April 26, 2017 Published: May 16, 2017 5603

DOI: 10.1021/acs.langmuir.7b00811 Langmuir 2017, 33, 5603−5612

Article

Langmuir tension fluids.28 Similarly, few numerical methodologies have been implemented to model the geometrical shape of pendant droplets and predict the critical Bo for droplet detachment from the surface.29−32 Lubarda and Talke33 determined the shape of a uniformly rotating droplet in gravity by iterative numerical integration of the governing nonlinear differential equations via variational analysis. Others have also worked on shape analysis of uniformly rotating droplets34−36 and free surface shape determination of liquid in a uniformly rotating cylinder.37 Related to shape prediction is the problem of shape analysis for determining interfacial properties such as surface tension and contact angle, among others. Many studies have been performed to measure the liquid−vapor interfacial tension,38 of which the pendant drop method39−43 and the sessile drop method44,45 are widely used. These methods rely on the extraction of the droplet profile from an experimental image and iteratively solve the Young−Laplace equation and estimate the physical parameters that describe the extracted drop profile most closely. With the advances in computational methods, droplet shape analysis techniques utilizing digital image processing46 such as axisymmetric drop shape analysis (ADSA)47 have been widely used48−50 over the past few decades. In this work, we present a comprehensive axisymmetric representation of the liquid−vapor interface using a vector parametrized cubic spline. The third-order spline representation used here is required to accurately capture the variations in curvature along the interface of a droplet. A thermodynamic free energy minimization based heuristic is then implemented to guide the iterative search and solve the Young−Laplace equation.5,51 We benchmark our approach by predicting the shapes of sessile and pendant droplets under the action of gravity. We show that the approach is also able to predict the shape of rotating droplets under the action of centrifugal force. At the end, the spline based representation is used to perform the inverse analysis, i.e., analyze a given droplet to predict the interfacial properties. We believe that simple and computationally inexpensive scheme proposed in this work has the potential to provide valuable insights into the physics of wetting required to enable fine fluidic manipulations in various applications. In addition, the simple inverse calculation scheme can serve as a powerful tool to improve applications such as surface tensiometry.

Figure 1. Schematic representation of the evolution of an arbitrary axisymmetric droplet shape (far from equilibrium shape) to the final equilibrium spherical cap shape.

interface until the pressure inside is uniform, i.e., the wellknown spherical/circular cap interface shape (red solid curve) with a constant curvature (χmean) is attained. It should be noted that this evolution strategy (heuristic) essentially follows thermodynamic free energy minimization via surface area minimization in this case. The algorithm is computationally inexpensive solution to the minimal surface problem52 which traditionally involves cumbersome variational formulation. Spline Based Representation. In this section, we discuss the geometric framework adopted for the representation of the liquid−vapor interface of droplets in this study. The interface is represented as a vector parametrized cubic spline S = [x(u) y(u)], where u ∈ [0, 1] represents the normalized spline parameter. The third-order spline representation used here over the usual constant curvature circular/spherical cap representation is required to accurately capture the variation in the magnitude and curvature of static droplets under the action of various types of forces such as gravity (Bo ≠ 0) or centrifugal. The details on the implementation of a cubic spline for modeling droplet interfaces can be found in the Supporting Information section S1. The spline-based parametrized representation makes it highly convenient to compute the contact angle, curvature, surface area, and volume of an axisymmetric droplet profile. All computations are performed using the commercial software package MATLAB.53 In order to simulate the equilibrium droplet shape, coordinates of the key points on the boundary of initial shape are used to generate the normalized parameter u. The generated parameter u is used to implement cubic spline in MATLAB53 over the curve. The left and the right contact angles are estimated by computing the inverse tangent of dy(u)/dx(u)|u=0,1, where u of 0 and 1 represent the start and the end points (contact points) of the parametrized spline. The principal curvatures at all points along the liquid−vapor interface are accordingly calculated as follows:

■

MODELING Our thermodynamic free energy minimization algorithm relies on the observation that equilibrium implies no liquid movement/circulation inside the droplet. In order to explain this, we assume the case of an arbitrarily axisymmetric pinned droplet (black solid curve) in the absence of any external force (Bo = 0) as shown in Figure 1. It is well-known that this shape does not reflect the equilibrium shape, which is that of a spherical cap (red solid curve). For example, the internal pressure PA on the liquid side next to the point A (blue square) is positive and equal to −2σχA,mean (we define mean geometric curvature χA,mean as −ve, concave outward) while the internal pressure PB on the liquid side next to the point B (red circle) is equal to −2σχB,mean, i.e., negative (the mean geometric curvature χB,mean is defined as +ve, concave inward). As a result, the natural tendency is to allow fluid flow from A to B, resulting in some intermediate shape schematically represented by the blue dashed curve. The droplet in principle will continue to allow such liquid redistribution (Video S1) all along the

χ1 =

x′y″ − y′x″ 1 = r1 (x′2 + y′2 )3/2

(1)

y′ 1 = 2 r2 x(x′ + y′2 )1/2

(2)

and χ2 =

where all primes represent the derivative with respect to the parameter u. Here r1 and r2 are the corresponding principal 5604


Article

Langmuir

Figure 2. Evolution of an arbitrarily shaped droplet (initial) with Bo = 0 to the final equilibrium shape on a surface with (a) pinned contact line. The evolution of corresponding χmean is shown in (b). Evolution of the same droplet to equilibrium, however, on a surface with (c) an advancing contact angle of 70°. The evolution of χmean for this case is shown in (d). The instantaneous value of local contact angles for each stage is shown in the legend.

relatively low energy intermediate shape (blue dotted curve) shown in Figure 1. Iterative use of these two operators eventually reduces the internal pressure to a constant, i.e., a spherical/circular cap droplet shape in equilibrium as discussed through Figure 1 (red solid curve). An example illustrating the implementation of this strategy (alternate application of PN and −PN operators) to attain equilibrium for an arbitrary shaped pinned droplet (Bo = 0) is presented in Figure 2a. The corresponding video is presented in the Supporting Information (Video S1). The blue squares and the red circles in Figure 2 and Video S1 represent the instantaneous locations of maximum and minimum curvature, i.e., the candidate perturbation points. The shape of the droplet is observed to evolve such that the variation in mean curvature χmean (pressure) is gradually minimized (intermediate shape), eventually reducing to a constant (brown dash-dot curve in Figure 2b). Please recall that u = 0 and u = 1 in Figure 2b represent the start and the end points (contact points), respectively, of the parametrized spline. We now explain the use of the DPN operator to simulate the evolution of the shape of advancing droplet on a surface with Young’s contact angle θ Y . The point on the spline corresponding to the contact line (u = 0 and u = 1) is incrementally moved outward (along the substrate) such that the instantaneous contact angle θ is reduced to a value infinitesimally smaller than the Young’s contact angle θY. However, such a movement also results in an artificial increase in the volume of the droplet (V → V + ΔV). Accordingly, this DPN operator is followed by a −PN operator such that the resulting volume increase is negated by the removal of ΔV from

radius of curvatures. The mean curvature (Figure 1) is defined as χmean =

χ1 + χ2 2

=

1⎛1 1⎞ ⎜ + ⎟ 2 ⎝ r1 r2 ⎠

(3)

It can be seen from eqs 1−3 that χmean is defined −ve for concave outward (point A in Figure 1) and +ve for concave inward (point B in Figure 1). The perimeter defined as Plv = ∫ u=1 u=0 dS is evaluated and is further used to compute the surface area and volume of the droplet using the Pappus’s centroid theorem.54 Perturbation Operators. Equipped with the spline based representation for liquid−vapor interface, we now define two unique perturbation operators, namely pinning (PN) and depinning (DPN), which when used iteratively will be shown to drive any initially assumed arbitrary axisymmetric droplet shape (for example, the black shape in Figure 1) to the final minimum energy equilibrium shape (for example, the red shape in Figure 1). As the name suggests, the contact line does not move during the PN operator. A small volume ΔV of liquid is added to the droplet profile at the location of the highest instantaneous mean curvature (minimum internal pressure, point B in Figure 1). If this operator is followed by −PN, i.e., a PN operator wherein same amount of liquid is removed (−ΔV) from the location of the lowest instantaneous mean curvature (maximum internal pressure, point A in Figure 1), then the overall volume of the droplet is maintained constant. The aforementioned two operators PN and −PN in turn mimic the natural tendency of the fluid to flow from the location of higher internal pressure (A) to lower pressure (B) and attain a 5605


Article

Langmuir

Figure 3. Methodology adopted for free energy minimization of the droplet.

shapes discussed in Video S1 and Video S2. The algorithm once again conveniently guides these randomly perturbed higher order shapes to the final equilibrium shape. A flowchart outlining the steps in the algorithm is shown in Figure 3. Please note that the simulation of the pinned droplet in Figure 2a essentially implies a very large advancing contact angle such that the condition “3” is never satisfied and the STICK loop (steps 2 → 3 → 4b → 5 → 2) is continuously invoked until equilibrium is attained. Moreover, if the initial nonspherical shape is defined such that the local contact angles are larger than the final contact angle of the shape in equilibrium (spherical cap), the contact line will recede in order to evolve to equilibrium. This scenario of an evolving droplet with receding contact line can be easily executed if in step 3 of Figure 3, “greater” and “advancing” are replaced by “smaller” and “receding”, respectively, and DPN and −PN in step 4a are replaced by −DPN and PN operators. Our algorithm is simple and provides an alternate solution to the Young−Laplace equation, bypassing the onerous task of solving partial differential equation. Discussion in Supporting Information section S2 presents a preliminary complexity analysis of our algorithm along with a discussion on the accuracy of algorithm. A detailed theoretical complexity analysis will be performed in the future.

the location of the lowest instantaneous mean curvature. These two operators (DPN followed by −PN) together mimic the natural tendency of the contact line to advance. The receding phase can similarly be mimicked by calling the −DPN (contact line is moved inward and volume of the droplet decreases artificially from V → V − ΔV) and PN (to negate the volume decrease) operators in turn. An example illustrating the use of PN and DPN operators to attain equilibrium for the arbitrarily shaped droplet on a surface with an advancing contact angle of 70° is shown in Figure 2c. The corresponding video (Video S2) is presented in the Supporting Information. Since the initial contact angle θ = 20° is lower than the advancing contact angle of 70°, the operators PN and −PN are called upon, and the droplet shape evolves toward equilibrium (intermediate I), similar to the example show in Figure 2a. The interface shape continues to evolve until the local contact angle becomes infinitesimally larger than the advancing contact angle of 70° (intermediate II). Now, the natural tendency of the contact line to advance is mimicked by invoking DPN operator along with the −PN operator. The outward movement of the contact line again lowers the contact angle. As a result, PN and −PN operators are called upon once again. This procedure wherein PN and −PN operators (STICK) are invoked if the local contact angle is incrementally lower than the advancing contact angle of the surface and DPN and −PN operators (SLIP) are invoked once the local contact angle exceeds the advancing value, if continued, eventually reduces the shape to the equilibrium spherical/circular cap (χmean = constant, brown dash-dot curve in Figure 2d) with a contact angle of 70° (red solid curve in Figure 2c and Video S2). Please note that our solution methodology is not limited up to fourth-order initial shapes as shown in Figure 2a. Since the splines used in our work are formed using piecewise polynomials, they enable us to simulate higher order shapes as well. The reader is referred to Video S3 and Video S4 wherein the order has been increased by randomly perturbing the initial

■

RESULTS

We now use the simulation strategy along with the cubic spline based interface representation to predict the equilibrium shapes of droplets which are influenced by gravity and centrifugal forces. Equilibrium Shape of Sessile Droplets (Bo ≠ 0). In this section, we used the perturbation operators to simulate the shape of droplets with various contact angles and different Bo. In contrast to the shapes simulated in Figure 2 where gravity effects were neglected and the liquid−vapor interface in 5606


Article

Langmuir

Figure 4. Evolution of an axisymmetric droplet shape (black, Bo = 0) with advancing contact angle (a) 60°, (b) 90°, (c) 120°, and (d) 150° to the final equilibrium shapes at various Bo.

Figure 5. Evolution of an axis symmetric droplet (black, Bo = 0) with receding contact angle (a) 60°, (b) 90°, (c) 120°, and (d) 150° to the final equilibrium shape at various Bo.

Example calculation for a sessile (Bo = 0.5, dashed curve) droplet for an advancing contact angle of 60° is shown in Figure 4a. The simulation starts with an initially assumed spherical cap interface shape shown using black solid curves. The sides of the sessile droplet bulge out due to gravity, the instantaneous local contact angle attains a value infinitesimally larger than 60°, DPN operator is invoked to advance the contact line, and the contact angle is reduced back to 60°.

equilibrium maintained a spherical cap shape, the variation of hydrostatic pressure (ρgh) inside large droplets cannot be neglected in comparison to the capillary pressure (2σχmean). Accordingly, fluid is redistributed such that the droplet attains a shape with a constant equivalent curvature (χeq = χmean − ρgh/ 2σ) and the droplet surface deviates from the spherical cap like geometry. Please note that the acceleration due to gravity g is considered positive in the downward direction. 5607


Article

Langmuir

Figure 6. Evolution of an axisymmetric droplet (black, Bo′ = 0) with an (a) advancing contact angle 56° to the final equilibrium shape at Bo′ = 5, Boω′ = 0 (nonrotating, blue dashed line) and Bo′ = 10, Boω′ = 15 (uniformly rotating, red dotted line). The volume V = 46.5 mm3 in both cases, while b = 6.10 mm and b = 8.63 mm, respectively. Evolution of a droplet with an (b) advancing contact angle 122° to the final equilibrium shape at Bo′ = 10, Boω′ = 0 (nonrotating, blue dashed line) and Bo′ = 20, Boω′ = 18.58 (uniformly rotating, red dotted line). The volume V = 275 mm3 in both cases, while b = 8.63 mm and b = 12.21 mm, respectively.

Figure 7. Evolution of a uniformly rotating axisymmetric droplet (black, Bo = 0, Boω = 0) with advancing contact angle (a) 60°, (b) 90°, (c) 120°, and (d) 150° to the final equilibrium shape at (a) Boω = 0.5, (b) Boω = 1, (c) Boω = 2, and (d) Boω = 3 and various Bo.

These operations minimize the variation in χeq along the liquid−vapor interface. Iterative use of these operations eventually reduces the droplet shape to equilibrium, as presented in Figure 4a. Results from similar analysis for droplets with different values of Bo on surfaces with advancing contact angles of 90°, 120°, and 150° are presented in Figures 4b−d, respectively. Respective plots of χmean and χeq for the results shown in Figure 4 are presented in Figure S4. The resulting shapes obtained via this exercise have been validated to be in excellent agreement with the results of Dimitrakopolous et al.19

A representative video for the advancing contact angle of 120°, and with a wide range of Bond numbers (Bo from 0 to 10 at an interval of 1), is presented in the Supporting Information (Video S5). Similar to Video S1 and Video S2, the blue squares and the red circles in the video represent the candidate perturbation points, i.e., the instantaneous locations of maximum and minimum equivalent curvatures χeq, respectively. Also, u = 0 and u = 1 represent the start and the end points (contact points) of the parametrized spline. Please note that while χeq becomes a constant throughout the liquid−vapor interface at equilibrium (refer to Supporting Information 5608


Article

Langmuir

Figure 8. (a) Shape of an axisymmetric droplet of a Bo = 5 on surface with advancing contact angle 60°. (b) σχeq is approaching minimum value at Bo = 5 (red triangle).

The simplicity of the proposed algorithm allows predicting the droplet shapes at various combinations of gravitational and centrifugal forces quite easily. For example, the equilibrium shapes of uniformly rotating droplet for different values of rotational Bond numbers (defined as Boω = ρω2L3/σ, where L is the characteristic dimension of the droplet) on surfaces with advancing contact angles of 60°, 90°, 120°, and 150° and in the absence of gravity are shown in Figure S6. Plots of χmean and χeq for the results shown in Figure S6 are presented in Figure S7. Video S7 shows the corresponding simulation on a surface with advancing contact angle of 150° with Boω values ranging from 0.5 to 3 at an interval of 0.5. Similarly, Figure 7 illustrates the equilibrium shape of a uniformly rotating pendant droplet (−ve gravity) with different values of rotational Bond numbers (Boω). The corresponding plots of χmean and χeq are shown in Figure S8 while a representative movie for the case of a pendant droplet with Boω = 1.5 and gravitational Bond number Bo ranging from −0.25 to −1 is presented through Video S8. Please note that while this surface does not have contact angle hysteresis wherein the advancing and receding contact angles are 90°, our simulation algorithm does not pose any limitations, and similar analysis on surfaces with contact angle hysteresis is also possible, but not shown here for brevity. Inverse Analysis. We now extend the application of spline based droplet shape representation and our thermodynamics based free energy minimization approach to develop a novel droplet shape analysis technique. This is first illustrated through the inverse or back calculation of Bond number Bo given the shape of liquid−vapor interface of an axisymmetric droplet. The inverse analysis relies on a simple strategy wherein the given profile is first digitized to obtain the coordinates of the points on the given liquid−vapor interface. As earlier, the obtained profile is then represented using spline and the mean curvature χmean is estimated at all the nodes. We then assume a value of Bo to estimate the equivalent curvature χeq = χmean − ρgh/2σ at all nodes, which is further used to estimate the standard deviation in the equivalent curvature (σχeq) for this assumed value of Bo. The whole process is repeated iteratively improving the guessed value of Bo which is allowed to span the expected parametric space of the droplet surface. σχeq is then plotted against the respective Bo to identify the Bo corresponding to the minimum value of σχeq, i.e., a fairly constant χeq implying the minimum energy equilibrium state. An example calculation for the droplet shape corresponding to Bo = 5 (data in Figure 4a) is shown in Figure 8. It can be seen that the standard deviation σχeq is minimized at Bo = 5, suggesting the ability of the inverse calculation to perform shape analysis as required in applications like tensiometry. It should be noted that for most cases we can identify the

section S3 for convergence criteria details), the geometric/ actual curvature χmean does not (Figure S4). Equilibrium Shape of Pendant Droplets (Bo ≠ 0). Similar to the sessile droplet simulations presented in the previous section, the simulation of pendant droplet is also started with an initially assumed spherical cap interface shape. Because of the −ve contribution of the hydrostatic pressure in this case, fluid is redistributed from the sides (the contact line would recede) to the apex in order to minimize the variation in χeq. We invoke −DPN and PN operators to recede the contact line such that sides of the droplet shrink inward and the instantaneous contact angle attains value larger than the receding contact angle. Further PN and −PN operators are invoked until instantaneous contact angle attains value incrementally smaller than the receding contact angle. This is again followed by invocation of −DPN and PN operators to recede the contact line. Iterative use of these operators eventually evolved the droplet shape to equilibrium. This is shown in Figure 5 for different values of Bond numbers Bo on surfaces with receding contact angles of 60°, 90°, 120°, and 150°. This is also apparent from Video S6 for receding contact angle of 90° with Bond number values ranging from −0.1 to −0.7 at an interval of −0.1. These shapes have also been validated to be in excellent agreement with the results of Dimitrakopolous et al.19 Plots of χmean and χeq for the results shown in Figure 5 are presented in Figure S5. Rotating Droplets. We next extend the energy minimization algorithm to simulate the equilibrium shapes of axisymmetric droplets rotating about the central axis passing through the apex, i.e., droplets subjected to centrifugal forces, where the directions of centrifugal and gravitational accelerations are perpendicular to each other. An additional pressure term ρω2x2/2 due to the centrifugal force is introduced in the ρgh

ρω2x 2

equation for equivalent curvature χeq = χmean − 2σ + 4σ , where ω is the angular velocity of the rotating droplet and x the horizontal distance of any point on the liquid vapor interface from the axis of rotation. Similar to droplets under the action of gravity, simulation is started with initial spherical cap interface shape, and the operators defined earlier are invoked iteratively until the equivalent mean curvature χeq is reduced to a constant all along the liquid−vapor interface of the droplet. The resulting shapes are presented in Figure 6 and have been validated to be in excellent agreement with the results of Lubarda and Talke.33 Table S3 includes coordinates of theoretical shape obtained via our algorithm and the digitized coordinates of droplet shape from the published result of Lubarda and Talke33 for Bo′ = 10 and Boω′ = 15. Lubarda and Talke33 defined gravitational and rotational Bond numbers as Bo′ = ρgb2/σ and Boω′ = ρω2b3/σ, respectively, where b is the radii of curvature at the apex of the droplet. 5609


Article

Langmuir

Figure 9. Inverse analysis of an axisymmetric droplet of a Bo = −0.2 and Boω = 3 on a surface with advancing contact angle 150°. (a) White hollow square corresponds minimum value of σχeq at Boω = 1. (b) Black solid square corresponds minimum value of σχeq at Boω = 4. (c) σmin χeq attains the minimum value at Boω = 3 (red triangle). (d) σχeq attains the minimum value at Bo = −0.2 corresponding to the minimum value of σmin χeq at Boω = 3 (red triangle).

in practical applications, and the variation of χeq can easily be checked to confirm the convergence of the solution algorithm. We believe that appropriate formulation for the parameter χeq can easily be developed to incorporate the effect of various other external forces including electric55,56 and acoustic.57 Moreover, the solution algorithm can be generalized to tackle any three-dimensional droplet shape either by using tensor product surfaces of cubic splines58 or a NURBS59 approximation based surface reconstruction algorithm.

physically realizable range of the guessed values of Bo by looking at the shape of the given droplet interface geometry, i.e., positive guess value of Bo for sessile droplet and negative guess value of Bo for pendant droplet. Since this method is computationally inexpensive, we run the algorithm over a wide range of Bo (e.g., say Bo ranging from −1000 to +1000) in the case of an ambiguity. The variation of σχeq with the guessed Bo is analyzed to narrow the range to a smaller zone until the instant of minimum error is identified accurately. For a uniformly rotating droplet under gravitational force, Boω is estimated along with Bo to estimate the equivalent ρgh

■

CONCLUSIONS We developed a vector valued parametrized cubic spline based representation for modeling liquid−vapor interface of axisymmetric droplets. The implementation mimics the liquid−vapor interface with a simplistic geometric optimization algorithm which equalizes the equivalent curvature and hence the pressure inside the droplet to minimize thermodynamic free energy and attain the equilibrium shape. The expression for equivalent curvature is fundamentally an augmented Young− Laplace equation which includes the effects of gravity and centrifugal forces in addition to the surface tension. We show that this strategy is simple and robust and is able to capture the shapes of sessile and pendant droplets over an extensive span of surface contact angles, and gravitational and rotation bond number values. The terms corresponding to surface tension, gravity, and centrifugal forces are independent of each other, making it convenient to perform the inverse analysis as well. We use droplet shapes from our analysis to show that we can accurately back-calculate gravitational and centrifugal Bond numbers independently. We believe that the proposed fluid− fluid interface modeling approach can tackle complicated droplet shapes encountered in various engineering applications, as desired.

ρω2x 2

curvature such that χeq = χmean − 2σ + 4σ is a constant at all nodes. We show this through the analysis of a droplet with Boω = 3 and Bo = −0.2 (dashed dotted curve in Figure 7d). We assume a value of Boω and estimate the standard deviation in equivalent curvature (σχeq) by changing the guess values of Bo. Sample results for Boω = 1 and Boω = 4 are shown in Figures 9a and 9b, respectively. The respective minimum values of σχeq are then identified (white hollow square for Figure 9a and black solid square for Figure 9b) and plotted with the value of Boω in Figure 9c. The minimum of σχeq (σmin χeq represented by red triangle in Figure 9c) corresponds to the correct solution, i.e., centrifugal, Boω = 3. Taking this established value of Boω, we revisit the plot of σχeq versus Bo to identify the Bo corresponding to the minimum value of σχeq (red triangle in Figure 9d) as the correct gravitational Bond number (Bo = −0.2). The minimum value of σχeq at various Bo and Boω for the data presented in Figure 9 is tabulated in Table S4. The examples illustrated above demonstrate that the droplet shape prediction and analysis technique developed in this work is simple, intuitive, and computationally inexpensive. Our approach works over a wide range of parameters encountered 5610


Article

Langmuir

■

(10) Rotenberg, Y.; Boruvka, L.; Neumann, A. W. The shape of nonaxisymmetric drops on inclined planar surfaces. J. Colloid Interface Sci. 1984, 102 (2), 424−434. (11) Feng, J. Q.; Basaran, O. A. Shear flow over a translationally symmetric cylindrical bubble pinned on a slot in a plane wall. J. Fluid Mech. 1994, 275, 351−378. (12) Muldowney, G. P.; Higdon, J. J. A spectral boundary element approach to three-dimensional Stokes flow. J. Fluid Mech. 1995, 298, 167−192. (13) Li, X.; Pozrikidis, C. Shear flow over a liquid drop adhering to a solid surface. J. Fluid Mech. 1996, 307, 167−190. (14) Dimitrakopoulos, P.; Higdon, J. J. Displacement of fluid droplets from solid surfaces in low-Reynolds-number shear flows. J. Fluid Mech. 1997, 336, 351−378. (15) Schleizer, A. D.; Bonnecaze, R. T. Displacement of a twodimensional immiscible droplet adhering to a wall in shear and pressure-driven flows. J. Fluid Mech. 1999, 383, 29−54. (16) Yon, S.; Pozrikidis, C. Deformation of a liquid drop adhering to a plane wall: Significance of the drop viscosity and the effect of an insoluble surfactant. Phys. Fluids 1999, 11 (6), 1297−1308. (17) Dimitrakopoulos, P.; Higdon, J. J. On the displacement of fluid bridges from solid surfaces in viscous pressure-driven flows. Phys. Fluids 2003, 15 (10), 3255−3258. (18) Spelt, P. D. Shear flow past two-dimensional droplets pinned or moving on an adhering channel wall at moderate Reynolds numbers: a numerical study. J. Fluid Mech. 2006, 561, 439−463. (19) Dimitrakopoulos, P. Gravitational effects on the deformation of a droplet adhering to a horizontal solid surface in shear flow. Phys. Fluids 2007, 19 (12), 122105. (20) Renken, F. P.; Subbarayan, G. NURBS-based solutions to inverse boundary problems in droplet shape prediction. Computer Methods in Appl. Mech. Eng. 2000, 190 (11), 1391−1406. (21) Stalder, A. F.; Kulik, G.; Sage, D.; Barbieri, L.; Hoffmann, P. A snake-based approach to accurate determination of both contact points and contact angles. Colloids Surf., A 2006, 286 (1), 92−103. (22) Lubarda, V. A.; Talke, K. A. Analysis of the equilibrium droplet shape based on an ellipsoidal droplet model. Langmuir 2011, 27 (17), 10705−10713. (23) Yu, Y. J.; Jung, H. Y.; Cho, H. G. A new water droplet model using metaball in the gravitational field. Computer & Graphics 1999, 23 (2), 213−222. (24) Malladi, R.; Sethian, J. A.; Vemuri, B. C. Shape modeling with front propagation: A level set approach. IEEE Trans. Pattern Anal. Mach. Intell. 1995, 17 (2), 158−175. (25) López, J.; Hernández, J.; Gómez, P.; Faura, F. A volume of fluid method based on multidimensional advection and spline interface reconstruction. J. Comput. Phys. 2004, 195 (2), 718−742. (26) Extrand, C. W.; Moon, S. I. Contact angles of liquid drops on super hydrophobic surfaces: understanding the role of flattening of drops by gravity. Langmuir 2010, 26 (22), 17090−17099. (27) Stalder, A. F.; Melchior, T.; Müller, M.; Sage, D.; Blu, T.; Unser, M. Low-bond axisymmetric drop shape analysis for surface tension and contact angle measurements of sessile drops. Colloids Surf., A 2010, 364 (1), 72−81. (28) Rienstra, S. W. The shape of a sessile drop for small and large surface tension. J. Eng. Math. 1990, 24 (3), 193−202. (29) Schulkes, R. M. The evolution and bifurcation of a pendant drop. J. Fluid Mech. 1994, 278, 83−100. (30) Chatterjee, J. Shape analysis based critical Eotvos numbers for buoyancy induced partial detachment of oil drops from hydrophilic surfaces. Adv. Colloid Interface Sci. 2002, 99 (2), 163−179. (31) Olgac, U.; Kayaalp, A. D.; Muradoglu, M. Buoyancy-driven motion and breakup of viscous drops in constricted capillaries. Int. J. Multiphase Flow 2006, 32 (9), 1055−1071. (32) Lamorgese, A.; Mauri, R. On the Buoyancy-Driven Detachment of a Wall-Bound Pendant Drop: Results of Phase-Field Simulations. Chem. Eng. Trans. 2015, 43, 1849−1854.

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.langmuir.7b00811. Section S1: discussion of implementation of the vector parametrized cubic splines; section S2: discussion about complexity and accuracy of algorithm; section S3: description of convergence criterion to establish the equilibrium shape of a droplet; section S4: comparison of χmean and χeq of evolving sessile and pendant droplets included in the text for different contact angles; section S5: comparison of coordinates of theoretical shape of uniformly rotating sessile droplet with published results; section S6: shapes of uniformly rotating droplets for different contact angles; section S7: comparison of χmean and χeq for uniformly rotating sessile and pendant droplets included in the text for different contact angles; section S8: information required to explain Figure 9; section S9: description of Videos S1−S8 (PDF) Video S1 (AVI) Video S2 (AVI) Video S3 (AVI) Video S4 (AVI) Video S5 (AVI) Video S6 (AVI) Video S7 (AVI) Video S8 (AVI)

■

AUTHOR INFORMATION

Corresponding Author

*Phone +91-612-302-8166; e-mail [email protected] (R.R.). ORCID

Rishi Raj: 0000-0002-3653-4288 Present Address

A.C.: Multi-Scale Multi-Physics Computation Lab, University of Texas at El Paso, El Paso, TX. Notes

The authors declare no competing financial interest.

■

REFERENCES

(1) Bear, J. Dynamics of Flow in Porous Media; Dover: New York, 1972. (2) Hynes, R. O. Integrins: versatility, modulation, and signaling in cell adhesion. Cell 1992, 69 (1), 11−25. (3) Stone, H. A.; Stroock, A. D.; Ajdari, A. Engineering flows in small devices: microfluidics toward a lab-on-a-chip. Annu. Rev. Fluid Mech. 2004, 36, 381−411. (4) Cristini, V.; Tan, Y. C. Theory and numerical simulation of droplet dynamics in complex flowsa review. Lab Chip 2004, 4 (4), 257−264. (5) Young, T. An essay on the cohesion of fluids. Philos. Trans. R. Soc. London 1805, 95, 65−87. (6) Bashforth, F.; Adams, J. C. An Attempt to Test the Theories of Capillary Action by Comparing the Theoretical and Measured Forms of Drops of Fluid; Cambridge University Press: Cambridge, 1883. (7) Larkin, B. K. Numerical solution of the equation of capillarity. J. Colloid Interface Sci. 1967, 23 (3), 305−312. (8) Brown, R. A.; Orr, F. M.; Scriven, L. E. Static drop on an inclined plate: analysis by the finite element method. J. Colloid Interface Sci. 1980, 73 (1), 76−87. (9) Milinazzo, F.; Shinbrot, M. A numerical study of a drop on a vertical wall. J. Colloid Interface Sci. 1988, 121 (1), 254−264. 5611


Article

Langmuir (33) Lubarda, V. A.; Talke, K. A. Configurational forces and shape of a sessile droplet on a rotating solid substrate. Theor. Appl. Mech. 2012, 39 (1), 27−54. (34) Angeles, J.; Thomson, M. Profile determination of a drop of liquid under surface tension, gravity and centrifugal forces. Comput. Mech. 1989, 4 (5), 329−344. (35) Bush, J. W.; Stone, H. A.; Bloxham, J. The motion of an inviscid drop in a bounded rotating fluid. Phys. Fluids A 1992, 4 (6), 1142− 1147. (36) Danov, K. D.; Dimova, S. N.; Ivanov, T. B.; Novev, J. K. Shape analysis of a rotating axisymmetric drop in gravitational field: Comparison of numerical schemes for real-time data processing. Colloids Surf., A 2016, 489, 75−85. (37) Lubarda, V. A. The shape of a liquid surface in a uniformly rotating cylinder in the presence of surface tension. Acta Mechanica 2013, 224 (7), 1365. (38) Neumann, A. W.; Good, R. J. In Techniques of Measuring Contact Angles; Good, R. J., Stromberg, R. R., Eds.; Plenum: New York, 1979; Vol. 11, pp 31−91. (39) Andreas, J. M.; Hauser, E. A.; Tucker, W. B. Boundary tension by pendant drops. J. Phys. Chem. 1937, 42 (8), 1001−1019. (40) Stauffer, C. E. The measurement of surface tension by the pendant drop technique. J. Phys. Chem. 1965, 69 (6), 1933−1938. (41) Heinig, P.; Steffen, P.; Wurlitzer, S.; Fischer, T. M. Twodimensional pendant droplet tensiometry in a Langmuir monolayer. Langmuir 2001, 17 (21), 6633−6637. (42) Drelich, J.; Fang, C.; White, C. L. Measurement of interfacial tension in fluid-fluid systems. Encycl. Surf. Colloid Sci. 2002, 3, 3158− 3163. (43) Berry, J. D.; Neeson, M. J.; Dagastine, R. R.; Chan, D. Y.; Tabor, R. F. Measurement of surface and interfacial tension using pendant drop tensiometry. J. Colloid Interface Sci. 2015, 454, 226−237. (44) Smolders, C. A.; Duyvis, E. M. Contact angles; wetting and dewetting of mercury: Part I. A critical examination of surface tension measurement by the sessile drop method. Recueil des travaux chimiques des Pays-Bas 1961, 80 (6), 635−649. (45) Butler, J. N.; Bloom, B. H. A curve-fitting method for calculating interfacial tension from the shape of a sessile drop. Surf. Sci. 1966, 4 (1), 1−7. (46) Emelyanenko, A. M.; Boinovich, L. B. The role of discretization in video image processing of sessile and pendant drop profiles. Colloids Surf., A 2001, 189 (1), 197−202. (47) Rotenberg, Y.; Boruvka, L.; Neumann, A. W. Determination of surface tension and contact angle from the shapes of axisymmetric fluid interfaces. J. Colloid Interface Sci. 1983, 93 (1), 169−183. (48) Cheng, P.; Li, D.; Boruvka, L.; Rotenberg, Y.; Neumann, A. W. Automation of axisymmetric drop shape analysis for measurements of interfacial tensions and contact angles. Colloids Surf. 1990, 43 (2), 151−167. (49) Cabezas, M. G.; Bateni, A.; Montanero, J. M.; Neumann, A. W. A new method of image processing in the analysis of axisymmetric drop shapes. Colloids Surf., A 2005, 255 (1), 193−200. (50) Marchuk, I. V.; Cheverda, V. V.; Strizhak, P. A.; Kabov, O. A. Determination of surface tension and contact angle by the axisymmetric bubble and droplet shape analysis. Thermophys. Aeromech. 2015, 22 (3), 297−303. (51) Laplace, P. S. Suppléments au Livre X. Traité de Méchanique Céleste; Gauthier-Villars: Paris, 1805; Vol. 4. (52) Lichnewsky, A.; Temam, R. Pseudosolutions of the timedependent minimal surface problem. J. Differ. Equations 1978, 30 (3), 340−364. (53) The MathWorks Inc. MATLAB, R2014b; The MathWorks Inc.: Natick, MA, 2014. (54) Goodman, A. W.; Goodman, G. Generalizations of the theorems of Pappus. Am. Math. Monthly 1969, 76 (4), 355−366. (55) Basaran, O. A.; Scriven, L. E. Axisymmetric shapes and stability of pendant and sessile drops in an electric field. J. Colloid Interface Sci. 1990, 140 (1), 10−30.

(56) Lac, E.; Homsy, G. M. Axisymmetric deformation and stability of a viscous drop in a steady electric field. J. Fluid Mech. 2007, 590, 239−264. (57) Yarin, A. L.; Pfaffenlehner, M.; Tropea, C. On the acoustic levitation of droplets. J. Fluid Mech. 1998, 356, 65−91. (58) Knott, G. D. Interpolating Cubic Splines; Springer Science & Business Media: Germany, 2012; pp 193−210. (59) Piegl, L.; Tiller, W. The NURBS Book; Springer Science & Business Media: Germany, 2012.

5612


Spline Based Shape Prediction and Analysis of ... - ACS Publications

Recommend Documents