Langevin Approach to Modeling of Small Levitating Ordered Droplet

Jun 27, 2018 - However, the interactions, which control the structure of the cluster, are still not well understood. ... repulsion forces (defined by ...
0 downloads 0 Views 3MB Size
Letter Cite This: J. Phys. Chem. Lett. 2018, 9, 3834−3838

pubs.acs.org/JPCL

Langevin Approach to Modeling of Small Levitating Ordered Droplet Clusters Nurken E. Aktaev,† Alexander A. Fedorets,† Edward Y. Bormashenko,‡ and Michael Nosonovsky*,†,§ †

University of Tumen, 6 Volodarskogo Street, Tumen 625003, Russia Department of Chemical Engineering, Biotechnology and Materials, Engineering Sciences Faculty, Ariel University, Ariel, Israel 40700 § Mechanical Engineering, University of WisconsinMilwaukee, 3200 North Cramer Street, Milwaukee, Wisconsin 53211, United States

Downloaded via AUSTRALIAN NATL UNIV on July 26, 2018 at 06:15:55 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.



S Supporting Information *

ABSTRACT: A self-assembled cluster of microdroplets levitating over a heated water surface is a fascinating phenomenon with potential applications for microreactors and for chemical and biological analysis of small volumes of liquids. Recently, we suggested a method to synthesize a cluster with an arbitrary number of small monodisperse droplets. However, the interactions, which control the structure of the cluster, are still not well understood. Here we propose a Langevin computational model considering the aerodynamic forces between the droplets and random diffusion-like fluctuations. Characteristic length and time scales and scaling relationships of interactions are discussed. The model shows excellent agreement with experimental observations for a small number of droplets.

M

extended periods of time. 17 We conjecture that the phenomenon, when combined with a spectrographic study of the droplets’ content, can be used for rapid biochemical in situ analysis. In a previous study, we demonstrated that clusters with an arbitrary small number of monodisperse droplets produce stable stationary configurations depending on the number of droplets.18 These unique configurations correspond to the minimum potential energy of the system, which is determined by the force balance. A local area of reduced pressure over a heated spot generates forces (defined by the effective potential W), which drive droplets toward the center of the heated area (the “main” force). When distances between the droplets are small and comparable with their diameters, short-range repulsion forces (defined by the potential U) act between the drops. The droplets in the cluster move, and the potential U is nonstationary. Note that, unlike for electrically charged particles, the potential field is not generated by the droplets themselves, but the interaction between the droplets is caused by local areas of high pressure when the droplets are close to each other.

icrodroplet ensembles have attracted the attention of scientists in the past decade in view of their microfluidics,1,2 biological,3 and chemical4,5 applications. It has been demonstrated that arrays of monodisperse microdroplets can revolutionize both the scale and speed of biological screening3 and enable control of chemical reactions on the scale of micro-, nano-, pico-, and femtoliters.6 Lab-on-chip, microelectromechanical systems (MEMS), and bio-MEMS devices exploit useful properties of droplets’ ensembles.7,8 In addition, ordered arrays of water droplets (capillary clusters) play a crucial role in the process of breath-figures self-assembly.9,10 Various experimental techniques were used to create and control capillary clusters, including the use of patterned and gradient surfaces and application of electromagnetic fields.11−13 A particularly interesting area is self-assembled levitating microdroplet clusters.14 Our group has demonstrated that the self-assembled monolayer of hexagonally ordered microdroplets (10−200 μm diameter) is condensed over a locally heated surface of water levitating at heights of the same order as the droplet diameter.15 The microdroplets do not coalesce due to complex aerodynamic forces between them in the ascending vapor-air jet. Typical water temperatures are 60−95 °C, although the phenomenon has been observed also at lower temperatures of about 20 °C.16 Changing the temperature and its gradient enables controlling the number of droplets, their density, and size. Furthermore, using infrared irradiation, it is possible to suppress droplet growth and stabilize them for © 2018 American Chemical Society

Received: May 31, 2018 Accepted: June 27, 2018 Published: June 27, 2018 3834

DOI: 10.1021/acs.jpclett.8b01693 J. Phys. Chem. Lett. 2018, 9, 3834−3838

Letter

The Journal of Physical Chemistry Letters For a small number of droplets on the cluster, N, the structure deviates from a hexagonal one and it depends on N. Here we develop a model of small cluster structure formation using the Langevin approach to include both aerodynamic and diffusive interactions and compare it with the experimental data to investigate the transition from a disordered structure to the hexagonal one. Forces that control the structure of the cluster have aerodynamic origin. A local area of reduced pressure over a heated spot generates forces (defined by the potential W) that drive droplets toward the center of the heated area. In other words, the droplet attains the area of reduced pressure on the surface {x, y}, which corresponds to minimum values of the potential W given by W (x , y ) =

k 2 (x + y 2 ) 2

(1)

where k is a phenomenological stiffness parameter of the potential. Although the stiffness can be determined experimentally from the frequency of free vibrations ω of a particle with mass m as k = ω2m, in this study we will consider a superposition of several potentials and therefore will use a different approach to find the value of k. To describe numerically the interaction between the droplets i and j, we introduce the phenomenological repulsion potential l 2 o o (yi − yj )2 | o o o o (xi − xj) = α expm − − } o o o o γx γy o o n ~

Figure 1. Effective potential for two droplets in the cluster located at x1 = x and x2 = −x. The effective “main” force, repulsion force, and their superposition are shown.

Ui[xi(t ), yi (t ), xj(t ), yj (t )]

(2)

Because the distance between the droplets is on the order of their diameter, |x1 − x2| ≈ d, using eq 3a one can estimate γx ≈ d2 = 10−8 m2. A previous study15 reported that the repulsion force is on the order of Frep = 10−9 N as using eqs 3a and 3b α α one can estimate Frep ≈ γ d ≈ d and finally α ≈ 10−13 J and k

where xn and yn are the coordinates of nth droplet. The phenomenological parameters γx and γy, which specify the fuzziness of the potential in the x- and y-directions, along with the parameter α specifying the height of the potential field’s peak, can be estimated from experimental data. Figure 1 shows the potentials W and U for a stable cluster of two droplets as a function of the distance between them. Note that the phenomenological repulsion potential is different from the one introduced for the viscous potential flow.19,20 The motion of the droplet in the combined potential field formed by the global potential W and local U can be approximated as a random walk of a Brownian particle in the phase space of the generalized momentum p{px,py} and coordinate ξ{x,y}. The characteristic dimensions of droplets are much smaller than the capillary length; thus, droplets are considered perfectly spherical. To model the interactions in a cluster, we must first estimate the characteristic time and length scales of these interactions. Let us consider the 1D motion of a cluster of two droplets of diameterd = 10−4 m. The stationary structure of the cluster implies that the combined potential is at its minimum, which yields the force balance equation 2| l o α o o (x1 − x 2) o kx1 = (x1 − x 2) expm − } o o o o γx γ x n ~

| l (x − x )2 o o α o 2 o kx 2 = (x 2 − x1) expm − 1 } o o o o γx γx n ~

x

≈ 10−7 kg/s2. These estimations by the order of magnitude can be used to investigate the stability of the model with respect to the change of the droplets’ size and mass. The lifespan of the cluster, τ, is on the order of dozens of seconds, and by stabilizing the cluster, it can be expanded to dozens of minutes and hours.17 The characteristic time scale of the momentum relaxation depends on dissipative properties of the system, and it can be estimated as the ratio of the mass of the droplet to the dissipation coefficient of the vapor-air jet m ≈ 0.03 s τp ≈ 3πdμ (4) where μ = 2 × 10−5 Pa·s is the coefficient of dynamic viscosity of wet air. The characteristic time scale of the coordinate relaxation can be estimated as a ratio of the dissipation coefficient of the vapor-air jet to the stiffness τx ≈

3πdμ ≈ 0.2 s k

(5)

Consequently, the cluster has the following hierarchy of time scales

(3a)

τ ≫ τx ≫ τp

(6)

The characteristic time of the thermal relaxation is estimated

(3b)

as 3835

DOI: 10.1021/acs.jpclett.8b01693 J. Phys. Chem. Lett. 2018, 9, 3834−3838

Letter

The Journal of Physical Chemistry Letters

Figure 2. Numerical modeling (potential maps) and experimental results18 for small clusters.

τtherm ≈

d2 αw

Considering the stochastic effect is required to exclude quasistable configurations (for example, when all particles are at the same straight line). The results show qualitative comparison of the numerical modeling results (potential maps with different colors corresponding to energy values) with the experimental observations (images of the small clusters). The model exactly reproduces the experimental results for clusters with less than 15 droplets (Figure 2, potential maps for small clusters with the number of droplets from 1 to 16 are presented in the Supporting Information, Figure S1). One particular experimental observation is that the structure of the cluster does not depend on the size of the droplet (note that the droplet’s diameter was always d < 100 μm, which is much smaller than the capillary length; therefore, the gravity effect was insignificant). This is warranted if the phenomenological parameters γx and γy are scaled as γx ∝ γy ∝ d2. To describe the cluster’s structure and compare numerical modeling and experimental results, we use the notation introduced earlier in ref 18; thus, N is the number of droplets in the cluster, G is the central group of droplets, and S = G + L2 is the nucleus including the central group and adjacent (second) droplet layer L. For larger clusters with the number of droplets from 16 to 28, the numerically obtained structures deviated from the experimental data (although there is exact fit for N18 and N23). Figure 3 shows the dependency of the cluster nucleus type on the number of droplets in the cluster. We attribute this deviation for clusters with more than 15 droplets to the idealized axially symmetric shape of the repulsion potential U assumed during the numerical modeling. We observed that the axially symmetric approximation remains valid until the third external layer of droplets emerges. However, even a minor deviation of the potential shape from the axial symmetry affects

(7) −7

−1

where αw ≈ 1.5 × 10 m s is the thermal diffusivity of water, which yields τtherm ≈ 10−3 s, smaller than τx and τp. When collective (macroscopic) parameters of the system, such as phenomenological k, γx, γy, and α, are changing slowly in comparison with the microscopic variables, namely, droplet coordinates and momenta, the Langevin equation can be applied. Thus, the motion of a cluster of N droplets is described by the stochastic Langevin equations21−23 i l o ÄÅ ÉÑ jjjj o o ÅÅ dp ÑÑ jj o W+ ÅÅ j ÑÑ jj−∇m ÅÅ ÑÑ = jj o o ÅÅ ÑÑ jj o o o ÅÅÅÇ dξj ÑÑÑÖ jjj n jj j 0 k

2

| zyz ÅÄ o o zz Å dt ÑÉÑ Ä o É o ∑ Ui}oo −6πrμzzzzz ÅÅÅÅÅ p ÑÑÑÑÑ ÅÅÅÅÅ ψ (t )dt ÑÑÑÑÑ o Ñ zz × ÅÅÅ j ÑÑÑ + ÅÅ i=0 o o zz ÅÅ dt ÑÑ ÅÅÇ 0 ÑÑÑÖ i≠j ~ zz ÅÅ m ÑÑ zz Ç Ö 1/m z{ N

(8) −5

where j is the droplet’s number, μ = 2 × 10 Pa·s is the coefficient of dynamic viscosity of wet air, and m = 5 × 10−10 kg is the droplet’s mass. Equation 8 yielded a stable solution, predicting the formation of a quasi-stationary cluster, appearing after a finite simulation time. The stochastic interaction caused by various fluctuations including those in the environment, heat source, etc. is given by the white noise ⟨ψ (t )⟩ = 0 ⟨ψx(t )ψy(t ′)⟩ = 2Dpδxyδ(t − t ′)

(9) (10)

where the factor of 2 is due to the dispersion, Dp = 6πrμT is the diffusion coefficient in the momentum space, temperature T = 343 K is the typical water surface temperature,14−17 δxy is the Kronecker delta, and δ(t) is Dirac’s delta function. 3836

DOI: 10.1021/acs.jpclett.8b01693 J. Phys. Chem. Lett. 2018, 9, 3834−3838

Letter

The Journal of Physical Chemistry Letters

while newly created small droplets were flown away by the strong vapor-air jet. Consequently, the cluster contained only those droplets, which were generated at the initial stage. The numerical solution of eq 9 was performed with the Euler−Maruyama method using a special program written in the programming language C(compiler MinGW). The differentiation step was decreased down to the value at which the results became insensitive to the step size. The numerical value of the simulation step was fixed for all calculations within the framework of our study and equal to τ = 5 × 10−4 s. At the initial moment of time, all particles were located at the x-axis, with such distances between them that the repulsion forces were insignificant in comparison with the conservative force. At every modeling step, the change of the coordinates of every particle was calculated considering the change in the repulsion forces. The computational process was run until the change of particles’ coordinates became insignificant and fluctuations were unable to change the structure of the cluster. The fluctuations were modeled by the function ψ(t) = (2Dp)1/2 dwt, where dwt is the independent normally distributed stochastic differential with zero mathematical expectation and variance dt. In discrete form, the function ψ(t) is written as b(Dpτ), where b is a random number subjected to the Gaussian distribution with variance 2. The linear congruential generator (LCR Lehmer) was used to generate random numbers.

Figure 3. Cluster types in modeling and experimental results.

the potential map significantly. Thus, for the cluster N12 (Figure 4), an introduction of asymmetry changes the central



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpclett.8b01693.



Figure 4. Potential maps of the cluster N12 for symmetric and asymmetric potentials U.

Plots with calculated potentials for clusters from 1 to 16 droplets (PDF)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Tel.: +1-414-229-2816. ORCID

Alexander A. Fedorets: 0000-0001-6595-3927 Edward Y. Bormashenko: 0000-0003-1356-2486 Michael Nosonovsky: 0000-0003-0980-3670

group of the cluster from G3 to G4. However, no clusters with the nucleus G4S12 were observed in experiment.18 For clusters larger than N15, higher terms of potential approximations (rather than the second-order terms used in eqs 1 and 2) may play a role. More accurate approximations can be obtained by solving the fluid dynamics problem with a larger number of droplets. We expect that more realistic asymmetric potentials will allow prediction of the structure of all small clusters.

Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS The authors are grateful to the Ministry of Education and Science of the Russian Federation (Project No 3.8191.2017/ BCh) for partial financial support of the present study. N.EA. thanks the Russian Foundation for Basic Research (Project No 18-38-00232 mol_a) for financial support.



METHODS The experimental methods of small cluster generation were described in detail earlier.18 An experimental setup included a cuvette with a 200 ± 5 μm thick layer of distilled water heated locally by laser irradiation from the bottom to reach the water surface temperature of 50−70 °C. To generate clusters with a constant number of almost identical monodisperse droplets, the process was divided into two stages. At the initial assembly stage, the power output of the laser was maintained at a small level, which allowed the droplets to migrate toward the center of the heated spot. After the required number of droplets was obtained, the power output was abruptly increased, leading to a higher speed of the ascending vapor-air jet. The droplets generated at the initial stage were large enough to levitate,



REFERENCES

(1) Beatus, T.; Bar-Ziv, R. H.; Tlusty, T. The Physics of 2D Microfluidic Droplet Ensembles. Phys. Rep. 2012, 516, 103−145. (2) 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. (3) Agresti, J.; Antipov, E.; Abate, E.; Ahn, K.; Rowat, A.; Baret, J.; Marquez, M.; Klibanov, A.; Griffiths, A.; Weitz, D. Ultrahighthroughput Screening in Drop-based Microfluidics for Directed Evolution. Proc. Natl. Acad. Sci. U. S. A. 2010, 107, 4004−4009.

3837

DOI: 10.1021/acs.jpclett.8b01693 J. Phys. Chem. Lett. 2018, 9, 3834−3838

Letter

The Journal of Physical Chemistry Letters (4) Song, H.; Tice, J.; Ismagilov, R. A Microfluidic System for Controlling Reaction Networks in Time. Angew. Chem. 2003, 115, 792−796. (5) Taniguchi, T.; Torii, T.; Higuchi, T. Chemical Reactions in Microdroplets by Electrostatic Manipulation of Droplets in Liquid Media. Lab Chip 2002, 2, 19−23. (6) Basova, E.; Foret, F. Droplet Microfluidics in (Bio)chemical Analysis. Analyst 2015, 140, 22−38. (7) Joensson, H. N.; Samuels, M. L.; Brouzes, E. R.; Medkova, M.; Uhlén, M.; Link, D. R.; Andersson-Svahn, H. Detection and Analysis of Low-Abundance Cell-Surface Biomarkers Using Enzymatic Amplification in Microfluidic Droplets. Angew. Chem., Int. Ed. 2009, 48, 2518−2521. (8) Juul, S.; Harmsen, C.; Nielsen, M. J.; Stougaard, M.; Leong, R. W.; Knudsen, B. R.; Ho, Y.-P. Single Cell Enzyme Diagnosis on the Chip. In Proc. 2013 IEEE 26th International Conference on Micro Electro Mechanical Systems (MEMS-2013), January 20−24, 2013; IEEE: Taipei, 2013; pp 911−914. DOI: DOI: 10.1109/ MEMSYS.2013.6474392. (9) Munoz-Bonilla, A.; Fernán dez-García, M.; RodriguezHernandez, J. Towards Hierarchically Ordered Functional Porous Polymeric Surfaces Prepared by the Breath Figures Approach. Prog. Polym. Sci. 2014, 39, 510−554. (10) Bormashenko, E. Breath-figure Self-assembly, a Versatile Method of Manufacturing Membranes and Porous Structures: Physical, Chemical and Technological Aspects. Membranes 2017, 7, 45. (11) Zheng, B.; Tice, J. D.; Ismagilov, R. F. Formation of Arrayed Droplets by Soft Lithography and Two-phase Fluid Flow, and Application in Protein Crystallization. Adv. Mater. 2004, 16, 1365− 1368. (12) Lyuksyutov, I. F.; Naugle, D. G.; Rathnayaka, K. D.D. On-chip Manipulation of Levitated Femtodroplets. Appl. Phys. Lett. 2004, 85, 1817. (13) Pompano, R. R.; Liu, W.; Du, W.; Ismagilov, R. F. Microfluidics Using Spatially Defined Arrays of Droplets in One, Two, and Three Dimensions. Annu. Rev. Anal. Chem. 2011, 4, 59−81. (14) Fedorets, A. A. Droplet Cluster. JETP Lett. 2004, 79, 372−374. (15) Fedorets, A. A.; Frenkel, M.; Shulzinger, E.; Dombrovsky, L. A.; Bormashenko, E.; Nosonovsky, M. Self-assembled Levitating Clusters of Water Droplets: Pattern-formation and Stability. Sci. Rep. 2017, 7, 1888. (16) Fedorets, A. A.; Dombrovsky, L. A.; Ryumin, P. I. Expanding the Temperature Range for Generation of Droplet Clusters over the Locally Heated Water Surface. Int. J. Heat Mass Transfer 2017, 113, 1054−1058. (17) Dombrovsky, L. A.; Fedorets, A. A.; Medvedev, D. N. The use of Infrared Irradiation to Stabilize Levitating Clusters of Water Droplets. Infrared Phys. Technol. 2016, 75, 124−132. (18) Fedorets, A. A.; Frenkel, M.; Bormashenko, E.; Nosonovsky, M. Small Levitating Ordered Droplet Clusters: Stability, Symmetry, and Voronoi Entropy. J. Phys. Chem. Lett. 2017, 8, 5599−5602. (19) Joseph, D. D. Viscous Potential Flow. J. Fluid Mech. 2003, 479, 191−197. (20) Joseph, D. D. Potential Flow of Viscous Fluids: Historical Notes. Int. J. Multiphase Flow 2006, 32, 285−310. (21) Yan, C.-C. S.; Chepyala, S. R.; Yen, C.-M.; Hsu, C.-P. Efficient and Flexible Implementation of Langevin Simulation for Gene Burst Production. Sci. Rep. 2017, 7, 16851. (22) Belyi, V. V. Thomson Scattering in Inhomogeneous Plasmas: the Role of the Fluctuation-Dissipation Theorem. Sci. Rep. 2018, 8, 7946. (23) Pavlova, E. G.; Aktaev, N. E.; Gontchar, I. I. Modified Kramers Formulas for the Decay Rate in Agreement With Dynamical Modeling. Phys. A 2012, 391, 6084−6100.

3838

DOI: 10.1021/acs.jpclett.8b01693 J. Phys. Chem. Lett. 2018, 9, 3834−3838