a many-body dissipative particle dynamics study

a School of Aerospace Engineering and Applied Mechanics, Tongji ... b School of Automotive Engineering, Shanghai University of Engineering Science, ...
0 downloads 0 Views 3MB Size
Subscriber access provided by READING UNIV

Article

Tuning drop motion by chemical chessboard-patterned surfaces: a many-body dissipative particle dynamics study Chensen Lin, Shuo Chen, lanlan xiao, and yang liu Langmuir, Just Accepted Manuscript • DOI: 10.1021/acs.langmuir.7b04162 • Publication Date (Web): 01 Feb 2018 Downloaded from http://pubs.acs.org on February 1, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Langmuir is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

Tuning drop motion by chemical chessboardpatterned surfaces: a many-body dissipative particle dynamics study Chensen Lin 1, Shuo Chena,*, Lanlan Xiaob, Yang Liuc a

School of Aerospace Engineering and Applied Mechanics, Tongji University, Shanghai 200092, China

b

School of Automotive Engineering, Shanghai University of Engineering Science, Shanghai 201620, China

c

Department of Mechanical Engineering, The Hong Kong Polytechnic University, Hung Hom, Kowloon, Hong Kong

*

Corresponding author. E-mail address: [email protected] (Shuo Chen).

ACS Paragon Plus Environment

1

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 28

Abstract

Controlling the motion of liquid drops on solid surface has broad technological implications. In this study, the many-body dissipative particle dynamics (MDPD) was employed to study the drop behaviors on chemical chessboard-patterned surfaces formed by square or triangular tiles. The scaling relationship of the model was established based on the surface tension, viscosity, and density of a real fluid, and an improved contact angle measurement technique was introduced to the MDPD system. For drops on a horizontal plane with different tile sizes, the equilibrium morphology was examined. The critical Bond number, i.e. the critical dimensionless force which is required to depin the drop, was found strongly affected by the size and shape of the tiles. Once the droplet begins to move, the tile pattern and size strongly affect the velocity fluctuation while weakly affect the average velocity. Interestingly, besides the common straight forward path, two more route patterns (zigzag and oblique) were observed by only tuning the tile angle, indicating that the advancing routes of the drop may vary according to the tile angle. To the author’s knowledge, this phenomenon has not been reported in the literature. This study provides a valuable tool to explore the possibility of passively control of drop’s motion by energy-free chemical heterogeneous surfaces, thus is helpful for engineers to design a surface that could manipulate the drop motion without external energy.

Keywords Many-body dissipative particle dynamics; Drop motion; Chemical heterogeneous surface

ACS Paragon Plus Environment

2

Page 3 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

1

Introduction

Controlling the motion of liquid drops on open surfaces is one of the major problems to be solved for the realization of reliable open microfluidic chips. The possibility of using chemically heterogeneous surfaces formed by suitable lyophobic regions to guide wetting drops along certain directions has recently attracted much attention both theoretically and experimentally 1. When drops are deposited on a surface functionalized with stripes of alternating wettability, they may assume elongated shapes, which are characterized by different contact angles in directions perpendicular and parallel to the stripes

2,3

. Studies on the drops sliding on striped surfaces also

reported an anisotropic behavior, i.e., drops slide more easily along the alternating stripes than across them 2,4, and periodic variations in the contact angles, possibly accompanied by fluctuations in the drop velocity, take place along this across direction 4. Changing the relative width of the stripes also has a strong influence on the morphology and the dynamics of drop 5. Besides the linear chemical pattern (strips), some researchers 6 designed and investigated various surfaces decorated with circular, square, or triangular domains (tiles) arranged in different configurations. It was reported that stick-slip motion was observed on all pattered surfaces, although it is less pronounced than on striped surfaces. In the dynamic regime, the slope of Ca/Bo mainly depends on the static equilibrium contact angle and weakly on the actual surface pattern. Through experiments and simulation, it was also found that a slower motion often occurs on a square chessboard and stripes pattern, compared with that on a triangular chessboard and aligned squares with a certain Bond number. Compared with stripes, the drop motion on these complex surfaces is affected not only by the wettability, fractions of hydrophobic materials, but also by the tile shapes and tile arrangement. The chessboard patterned surface is one of the two dimensional periodic patterns with equal hydrophilic and hydrophobic domains. Better understanding the drop

ACS Paragon Plus Environment

3

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 28

behavior on chessboard surface lays a good foundation for other chemical heterogeneous surface. In this study, two variant chessboard models: square chessboard and triangular chessboard will be investigated. Simulation could offer great flexibility in investigating a variety of load conditions and performing local measurements of capillary, viscous, and body forces, which are otherwise difficult to obtain by experimental means 5. To simulate the sliding droplets on these surfaces, different approaches might be followed 7–11. Assuming the fluids as a continuum, one can solve the Navier-Stokes equation describing the evolution of the full velocity field complemented with appropriate boundary conditions at the solid-liquid and the liquid-air interfaces 10. Alternatively, one may combine the bulk equations with an additional phase-field dynamics, thereby modelling the liquid-air interface 12,13. More fundamental approaches are particle-based discrete methods, such as molecular dynamics (MD) 14 or dissipative particle dynamics (DPD) 15. These methods are capable of simulating phenomena where the continuum assumption breaks down. However, MD requires a prohibitively large computational effort and thus are restricted to nano-droplets 10. The advantage of dissipative particle dynamics (DPD) resides in the coarse-graining notion with particle interaction under a soft repulsive potential, making it a mesoscopic model of molecular dynamics 16. Therefore, a more sophisticated way to incorporate the essential microscopic or mesoscopic physics while recovering the macroscopic laws and properties at affordable computational cost. In single-species fluid problems, however, the standard DPD method has a fundamental limitation, in that the repulsive soft potential alone cannot reproduce surface tension of liquid/vapor interface. This potential leads to a predominantly quadratic pressure-density equation of state (EOS) 17, while a higher-order pressure-density curve is necessary for the coexistence of the liquid and vapor phases . The many-body DPD (MDPD) method includes an

ACS Paragon Plus Environment

4

Page 5 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

attractive force with the amplitude of the soft repulsion proportional to the local density of the particles, thus achieves a cubic pressure-density relation 18. The MDPD was extensively investigated and has been employed in studies of the contact angle in capillary flow 19, liquid droplets on solid surfaces 20, and bubble formation 21, etc. The paper is organized as follows: Section 2 deals with the numerical techniques and parameters applied in this study, the corresponding results are described in Section 3, while conclusions are reported in Section 4. All the simulations are enabled by the particle dynamics software code LAMMPS 22 with the addition of a homemade MDPD class.

2

Methods and conditions

2.1

MDPD scheme

The MDPD inherits the three pairwise-additive inter-particle forces formulation of the standard 𝐷𝐷 𝑅𝑅 + 𝐅𝐅𝑖𝑖𝑖𝑖 . The conservative, dissipative, and random forces of this DPD scheme, F𝑖𝑖 = ∑𝑖𝑖≠𝑗𝑗 𝐅𝐅𝑖𝑖𝑖𝑖𝐶𝐶 + 𝐅𝐅𝑖𝑖𝑖𝑖

expression are defined, respectively, as 𝐅𝐅𝑖𝑖𝑖𝑖𝐶𝐶 = 𝐹𝐹𝑖𝑖𝑖𝑖𝐶𝐶 (𝑟𝑟𝑖𝑖𝑖𝑖 )𝐫𝐫�𝑖𝑖𝑖𝑖

(1)

𝐷𝐷 𝐅𝐅𝑖𝑖𝑖𝑖 = −𝛾𝛾𝜔𝜔𝐷𝐷 (𝑟𝑟𝑖𝑖𝑖𝑖 )(𝐯𝐯𝑖𝑖𝑖𝑖 ⋅ 𝐫𝐫�𝑖𝑖𝑖𝑖 )𝐫𝐫�𝑖𝑖𝑖𝑖

(2)

𝑅𝑅 𝐅𝐅𝑖𝑖𝑖𝑖 = 𝜉𝜉𝜃𝜃𝑖𝑖𝑖𝑖 𝜔𝜔𝑅𝑅 ∆𝑡𝑡 −1/2 (𝑟𝑟𝑖𝑖𝑖𝑖 )𝐫𝐫�𝑖𝑖𝑖𝑖

(3)

where 𝐫𝐫�𝑖𝑖𝑖𝑖 = 𝐫𝐫𝑖𝑖𝑖𝑖 /𝑟𝑟𝑖𝑖𝑖𝑖 , 𝐯𝐯𝑖𝑖𝑖𝑖 = 𝐯𝐯𝑖𝑖 − 𝐯𝐯𝑗𝑗 and ∆𝑡𝑡 is the simulation time step.

Warren’s approach 23 is pursued for 𝐅𝐅𝑖𝑖𝑖𝑖𝐶𝐶 ; the repulsive force depends on a weighted average of the local density, whereas the attractive force is density-independent,

ACS Paragon Plus Environment

5

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

𝐹𝐹𝑖𝑖𝑖𝑖𝐶𝐶 = 𝐴𝐴𝑖𝑖𝑖𝑖 𝑤𝑤𝑐𝑐 (𝑟𝑟𝑖𝑖𝑖𝑖 ) + 𝐵𝐵𝑖𝑖𝑖𝑖 (𝜌𝜌̅𝑖𝑖 + 𝜌𝜌̅𝑗𝑗 )𝑤𝑤𝑑𝑑 (𝑟𝑟𝑖𝑖𝑖𝑖 )

Page 6 of 28

(4)

The weight functions 𝑤𝑤𝑐𝑐 (𝑟𝑟) = (1 − 𝑟𝑟/𝑟𝑟𝑐𝑐 ) and 𝑤𝑤𝑑𝑑 (𝑟𝑟) = (1 − 𝑟𝑟/𝑟𝑟𝑑𝑑 ) vanish for 𝑟𝑟 > 𝑟𝑟𝑐𝑐 and 𝑟𝑟 > 𝑟𝑟𝑑𝑑 ,

respectively. Since a DPD method with a single range may not have a stable interface 24, in Eq.(4) the repulsive contribution is set to act at a shorter range 𝑟𝑟𝑑𝑑 < 𝑟𝑟𝑐𝑐 than the soft pair attractive potential.

The many-body repulsion is chosen in the form of a self-energy per particle which is quadratic in the local density, 𝐵𝐵𝑖𝑖𝑖𝑖 (𝜌𝜌̅𝑖𝑖 + 𝜌𝜌̅𝑗𝑗 )𝑤𝑤𝑑𝑑 (𝑟𝑟𝑖𝑖𝑖𝑖 ), where 𝐵𝐵 > 0. The density for each particle is defined as 𝜌𝜌̅𝑖𝑖 = ∑𝑗𝑗≠𝑖𝑖 𝑤𝑤𝜌𝜌 (𝑟𝑟𝑖𝑖𝑖𝑖 )

(5)

and its weight function 𝑤𝑤𝜌𝜌 is defined as 𝑤𝑤𝜌𝜌 (𝑟𝑟) =

15

2𝜋𝜋𝑟𝑟𝑑𝑑3

(1 − 𝑟𝑟/𝑟𝑟𝑑𝑑 )2

(6)

𝑤𝑤𝜌𝜌 vanishes for 𝑟𝑟 > 𝑟𝑟𝑑𝑑 and for convenience is normalized.

The DPD thermostat consists of random and dissipative forces. The 𝜃𝜃𝑖𝑖𝑖𝑖 from Eq. (3) are

independent identically distributed Gaussian random numbers with zero mean and unit variance. The equilibrium temperature 𝑇𝑇 is maintained through the condition posed by the fluctuation-

dissipation theorem 𝜔𝜔𝐷𝐷 (𝑟𝑟) = [𝜔𝜔𝑅𝑅 (𝑟𝑟)]2

(7)

𝜉𝜉 2 = 2𝛾𝛾𝑘𝑘𝐵𝐵 𝑇𝑇

(8)

𝜔𝜔𝐷𝐷 (𝑟𝑟) = (1 − 𝑟𝑟/𝑟𝑟𝑐𝑐 )2

(9)

where 𝑘𝑘𝐵𝐵 is the Boltzmann constant. The weight function for the dissipative force is

ACS Paragon Plus Environment

6

Page 7 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

Details about the choice of the random and dissipative coefficients 𝜉𝜉 and 𝛾𝛾 can be found in Ref. 25

. All the simulations presented in this work were carried out with the modified velocity Verlet

algorithm 25 using the value of 0.5 for the empirical parameter. As in any other method in computational fluid dynamics, the issue of boundary conditions has to be addressed. There are many studies focused on boundary conditions in DPD. Many different boundary treatment techniques have been introduced for different types of problems. In this work, we use the method provided by Li et al. 26, which is verified capable of handling complex wall condition. In this method, by introducing an indicator variable of boundary volume fraction for each fluid particle, the boundary of arbitrary-shape objects is detected on–the-fly for the moving fluid particles. With the help of the modified effective dissipative interaction, the no-slip condition is satisfied and the fluctuation of fluid density near the wall is well controlled. 2.2

Simulation procedures

In experiments, chemical patterned surfaces can be produced by photolithography 27. Through AFM measurements, the roughness of surface is shown less than 10 nm both on printed (octadecyltrichlorosilane) and bare regions (glass) 6. With the AFM in topography mode, it is not possible to distinguish the printed areas from the bare glass regions, indicating the thickness of the printed layer is negligible with respect to the native roughness of the whole surface. Accordingly, in our simulation the roughness of the surface was made homogenous. To present chemical hetrogeneous, only different wettablities was assigned to wall particles in different wettability regions. Sliding simulations were performed depicted in Fig. 1(a). The driving froce of the drop was provided by the down-plane component of the drop weight, i.e. 𝐹𝐹𝑔𝑔 = 𝜌𝜌𝜌𝜌𝜌𝜌sin𝛼𝛼, acting on a drop ACS Paragon Plus Environment

7

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 28

of volume 𝑉𝑉 and density 𝜌𝜌 sliding dwon an inclined plane tilted by an angle 𝛼𝛼. Two surface

patterns were investigated in this study: square chessboard (SQ) and triangular chessboard (TR). The size of the tile is represented by the row spacing (W).

Fig. 1 (a) Sliding drop on an inclined plane tilted by an angle 𝛼𝛼. The surfaces is functionalized

with squares or triangles of different wettability with row spacing W. Red tile (dark) corresponds to hydrophobic domains; grey (light) tile corresponds to hydrophilic domains. (b)-(e) Drops of same size on surfaces with varied tile shape and size.

ACS Paragon Plus Environment

8

Page 9 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

2.3

Parametrization and contact angle measurement

We simulate drops of a 80% w/w glycerol/water solution (density 𝜌𝜌 = 1.21 g cm−3 , viscosity 𝜂𝜂 =

52 cP , and surface tension 𝛾𝛾 = 65.3 mN m−1 at 𝑇𝑇 = 23 ℃ ), which has been used in an

experimental study of the drop motion on simple chemical heterogeneous surface — a chemical

step 11. Unfortunately, tuning force parameters in MDPD to represent a specific kind of fluid is still open to discuss. A brief discussion on the relationship between the parameters, coarse-graining level and the properties could be found in Ref. 28. In our study, the interaction parameters set are described by 𝐴𝐴ll = −40, 𝐵𝐵ll = 25, 𝑟𝑟𝑑𝑑 = 0.75, which are proposed in Ref. 28 and thereafter being representative values for solid-liquid-vapor coexisting system

19,29

. The number density and the

surface tension coefficient can be calculated with numerical tests. We find 𝜌𝜌 = 6.08, and 𝛾𝛾 = 7.28 with minimal error compared to Ref. 11 (𝜌𝜌 = 6.10, and 𝛾𝛾 = 7.30).

Despite the convenience of referring to dimensionless quantities due to the convention that DPD operating in reduced units, it is instructive to map the physical units of length, mass, and time to a specific set of properties of a real liquid. In the following sections, quantity with square bracket denotes the scaling relationship of the quantity. For instance, 𝑙𝑙physical = [𝑙𝑙] ∗ 𝑙𝑙MDPD . By

correlating the density, viscosity, and surface tension of our MDPD fluid and 80% w/w glycerol/water solution, we obtain [𝑙𝑙] = 6.2349 × 10−5 m , [𝑚𝑚] = 4.8219 × 10−11 kg , [𝑡𝑡] = 7.3322 × 10−5 s. Several key properties in this study are shown in Table 1. The unitless and

physical properties could convert to each other according to the established scaling relationship. For simplicity and consistency with previous DPD research work, the simulation results are presented in reduced units in the rest of the paper.

ACS Paragon Plus Environment

9

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 28

Table 1 Scaling of MDPD Description

DPD units

Physical units

𝜌𝜌

density

6.08

1.21 g cm−3

𝑉𝑉

volume

164636

40 µl

𝛾𝛾

surface

𝑔𝑔

gravity

𝑟𝑟

radius

𝜂𝜂

viscosity

tension

34

4.93 7.28

8.45 × 10−4

2.12 mm 52 cP

65.3 mN m−1 9.8 m s−2

Contact angle (CA) is of vital importance in describing the surface wettability. In experimental study, the contact angle could be measured using contact angle measuring device (i.e. Dyne Theta optical tensiometer, Dataphysics OCA20, etc.). However, there is no standard procedure in particle simulations to measure the contact angle. A simplified and commonly adopted approach uses 3 characteristic points (i.e. point a, b, and c in Fig. 2(a)) on the drop to determine a circle whose tangent indicates the contact angle. This method is easy-to-implement, however, not capable of handling advancing and receding contact angle. What’s more, the hypothesis of drop shape being perfect hemispherical fails when the gravity is non-zero and the surface is anisotropic. In this study, we proposed an improved method to measure contact angle in particle system. To start with, the particles composing the sessile drop must be identified. This step is necessary when measuring the CA every time since the “condensation” and “evaporation” is constantly ongoing and the property of vapor and liquid particles are essentially identical. Technically the identification is made possible by grouping particles into different clusters, in which particles are within the cutoff

ACS Paragon Plus Environment

10

Page 11 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

distance from other particles in the same cluster. The cluster with the majority of particles in the system is recognized as the sessile drop. The results are shown by coloring in Fig. 2, the vapor particles are marked yellow, while the sessile drop particles are marked blue. Next, the surface particles of the sessile drop are identified by evaluating the number of their neighbors. Then a sub collection of the surface on the direction of interests is isolated for measuring the CA. As a thumb of rule, the arc-like sub collection is in made 5% drop diameter wide and 10% drop diameter high from the wall. Two specific arcs are colored red in Fig. 2(b) to calculate CA on the left and right. Final, the arc particles are used to determine a best fitting circle by least-squares approximation. Noteworthy, the CA calculated from one snapshot of the system exhibits noticeable fluctuation due to the thermal fluctuation of the surface. To obtain an experimental consistent CA, the averaging of 100 steps is used as the output CA in this study. The advantage of this improved CA measure method is the capability to access CA in arbitrary direction, while being insensitive to the deviation of drop shape from the perfect hemisphere.

Fig. 2 (a) Conventional contact angle measurement with three points to circle fitting (without gravity); (b) Improved contact angle meacurement technique (with gravity).

ACS Paragon Plus Environment

11

Langmuir

The wettability of the surface to the liquid could be tuned by changing the solid-liquid interaction parameter

19,30,31

. The correlation between the force parameter and the CA is investigated by

varying the 𝐴𝐴sl from 4 to 36. It is evident in Fig. 3 that the CA is responsive to 𝐴𝐴sl in the middle,

but insensitive when 𝐴𝐴sl continue growing or decreasing. By fitting these measured data, we have the following relationship to estimate the CA for a given 𝐴𝐴sl .

(10)

𝜃𝜃s = 74.03 + 72.11 ∗ sin(𝜋𝜋 ∗ (14.72 − 𝐴𝐴sl )/36.49)

In the following cases, the relatively hydrophilic domains are represented by particles with 𝐴𝐴sl = −20 (CA~87°), and hydrophobic domains 𝐴𝐴sl = −10 (CA~132°). 160

Measured contact angle Fit curve of measured contact angle

140

Static Contact Angle θs (°)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 28

120 100 80 60 40 20 0 0

-4

-8

-12

-16

-20

-24

-28

-32

-36

Solid-Liquid Interaction Parameter Asl

Fig. 3 The relationship between solid-liquid interaction parameter and static contact angle. The size of droplet strongly affects the dynamic behavior on surface considering that the driving force and the resistance force did not grow with the same order. Usually, a 40 µL droplet is typical

for study of drop motion on heterogeneous surface in experiments 6,11. According to the established scaling relationship, 40 µL corresponds to a drop with diameter 34 DPD units containing a million MDPD particles. This system size overwhelmed out current computing capacity. To develop a

ACS Paragon Plus Environment

12

Page 13 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

viable solution and take advantage of the flexibility of computer modeling, we take a close look at the problem. By varying 𝛼𝛼 between 5° and 45°, the in-plane component of the body force can be described by 2

a Bond number Bo = 𝜌𝜌𝜌𝜌𝑉𝑉 3 sin𝛼𝛼�𝛾𝛾 ranging between 0.2 and 1.5. The typical velocities of the drop

steadily sliding on hydrophilic and hydrophobic domains range from 0.01 to 0.1. Accordingly, the maximum Weber and Capillary number of the drop with dimension 𝐿𝐿0 ~𝑉𝑉 1/3 were We = 𝜌𝜌𝑈𝑈 2 𝐿𝐿0 ⁄𝛾𝛾 ~0.225 and Ca = 𝜂𝜂𝜂𝜂/𝛾𝛾~0.06, respectively. The close to unity Bo and small We and

Ca indicate that drop motion was mainly governed by an interplay of gravity and interfacial forces. In order to reproduce the case, the key factor Bo should maintain equivalent in simulation. If we

halve the drop radius (new r=17), while magnifying the gravity (new g=0.00338), the Bo number will be the same in experiment and simulation.

3 3.1

Results and discussion Morphology

A drop deposited on a horizontal homogeneous surface assumes a hemispherical-like shape and the apparent contact angle is the same as measured in any direction. Nevertheless, on chemical heterogeneous surface, the drop’s morphology and the CA both twist to adapt the alternation of the surface wettability. Some studies minimize these diversities by making the droplet covering plenty of chemical tiles, so as to make it possible to explore the drop’s dynamic behavior on chemical homogeneous and heterogeneous surfaces while nearly avoiding the influence of the wetting area deformation and CA fluctuation. However, the contact line and the dynamics of the drop would be obviously different if the tile size increases. Other studies focus on the deformation and CA deviation 3,32,33 of the drop on chemical heterogeneous surface. To the best

ACS Paragon Plus Environment

13

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 28

of our knowledge, these studies are all on chemical striped surfaces. In this section, we examine the morphology of the wetting area as the chemical pattern size grows. The wetting area is defined as the contact area between the drop and the surface. Its contour fundamentally influences the drop shape. In experiments, a clear view of the wetting area is hardly reachable if the drop’s projection is broader than the wetting area. Fortunately, we could easily remove the upper part of the drop in simulation. In Fig. 4(a), the wetting area is characterized by the thin particle layer below 1 DPD length from the surface, marked blue; the main body of the drop is made half transparent. These three view illustrate the wetting area is different from what we see from up down. Fig. 4(b)(c) demonstrate the wetting area on horizontal SQ and TR surfaces with different tile sizes. The main body of the drop is made invisible. It is evident that the wetting area’s morphology is various. It is roughly a circle while W is small, confirming that by reducing W, the effect of heterogeneous surface on wetting area could be diminished. The wetting area begins to be constrained by the tile shape and size as W increases. Different morphologies have been observed on different surfaces, such as, rectangular (SQ W=8≈0.47r), capsule-shape (SQ W=20≈1.18r), triangular (TR W=16≈0.94r), trapezoid (TR W=18≈1.06r), etc. Compared to wetting area on stripe-patterned surfaces, drop on chessboard-patterned surface exhibit further more diversity and attractiveness.

ACS Paragon Plus Environment

14

Page 15 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

Fig. 4 (a) Illustration of the of wetting area (Front, top, and perspective view); (b) Top view of wetting area on different size of SQ tiles; (c) Top view of wetting area on different size of TR tiles. Moreover, displacement of drop is also detected in some cases. The dripping position, which is marked by black cross in Fig. 4(b)(c), is standing on the cross point of hydrophilic and hydrophobic

ACS Paragon Plus Environment

15

Langmuir

domains. In most cases the drop will rest centered in the dripping position. In other cases, i.e. SQ W=10, the equilibrium position deviated from the initial position a remarkable distance. This could be explained by the theory that the final equilibrium wetting area has a lower energy than spreading around the dripping point. Interestingly, in the cases with small W, slight offset was also recognized. We believe this offset is caused by a different mechanism. By repeating the same configurations over 20 times, we found the offset direction and distance are hard to predict. A possible reason for this is the initial disturbance when dropped down, companied with the little constrain on the drop free movement from the small W surfaces.

1200

90 80

1000

70 800

60 50

600

40 400

30 Area of contact (SQ) Area of contact (TR) Hydrophilic fraction of contact area (SQ) Hydrophilic fraction of contact area (TR)

200

20 10

Hydrophilic fraction of contact area (%)

100

Area of contact

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 28

0

0 0

5

10

15

20

25

W

Fig. 5 The area of wetting and hydrophilic fraction of wetting area on surfaces with different W. To better understand the equilibrium state on surfaces with different tile size, we measure the area of wetting and the hydrophilic fraction of wetting area. As shown in Fig. 5, the area of wetting fluctuated, especially when W is close to the drop radius (r=17). On the other hand, we found the hydrophilic fraction value experience a relative steady slope. The fraction stands at 50% on small

ACS Paragon Plus Environment

16

Page 17 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

W surfaces, and could reach 70% (SQ) and 66% (TR) when the tile size is much bigger than the drop size. It is natural to embrace the hydrophilic surface, and with larger tiles, the drop has more free space to adjust itself for higher hydrophilic fraction. 3.2

Bond effect

Sliding measurements have been performed for several inclinations on SQ/TR patterned surfaces with W=5 and 10, respectively. From the time evolution of the drop position, the mean velocity U has been determined for different inclination angle 𝛼𝛼 . When the driving force based on the component of gravity parallel to the plane exceeds a threshold, the drops start moving at a steady velocity which increases linearly with Bo. The results obtained are summarized in Fig. 6 in terms of the dimensionless capillary number (Ca) versus Bond number (Bo). We found that in the dynamic region, the slope of ∆Ca/∆Bo depends weakly on the actual surface pattern size and

shape. This finding is consistent with the results of previous studies 6,27. However, the divergence slightly grows as Bo increase. Overall, surfaces with W=5 seems more slippery than W=10 in both SQ and TR cases; and TR pattern is more slippery than SQ pattern. The difference of CA is a result of the different energy dissipation level caused by pinning-depinning transition of the contact line on each specific surface. In other words, for smaller CA cases, larger part of driving energy is stored in the larger periodic deformations. Aside from the cases above, another series (TR W=5) is also provided with hollow symbols in Fig. 6. These additional cases apply the original parameter set (r=34, g=0.000845) disused in Section 2.3. Insignificant difference between the cases with original set and modified set (r=17, g=0.00338) validate the effectiveness of the modified parameter set.

ACS Paragon Plus Environment

17

Langmuir

0.16 0.14 0.12

TR W=5 TR W=5 (verification) TR W=10 SQ W=5 SQ W=10

0.10

Ca

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 28

0.08 0.06 0.04 0.02 0.00 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2 1.3

Bo

Fig. 6 Non-dimensional velocity of the drop as a function of plane inclination. Ca = 𝜂𝜂𝜂𝜂/𝛾𝛾 is the 2

capillary number based on the drop velocity 𝑈𝑈, and Bo = 𝜌𝜌𝜌𝜌𝑉𝑉 3 sin𝛼𝛼�𝛾𝛾. 3.3

Pinning

Despite the slope of ∆Ca/∆Bo in Fig. 6, it is also noteworthy that the onset of drop motion is strongly affected by the surface. The pinning phenomenon on chemical heterogeneous surfaces

has attracted a lot of attentions. For stripe patterned surfaces, the studies mainly concentrated on the fraction of different wettability regions with the same periodicity 5,34. For chessboard patterned surfaces, the fraction is naturally 50% for both SQ and TR, however, several other factors influence the pinning phenomenon. It has been reported that the arrangement of the tile have strong impact on the onset of drop motion 6. In this section, we focused on the effect of the tile size. To this end, we test the surface with W ranging from 2.5 to 40 to find out the critical Bo number (Boc ). The

results are plotted in Fig. 7. The curves exhibit obvious non-linear relationship between the tile

ACS Paragon Plus Environment

18

Page 19 of 28

size and the Boc . From the chart, we could tell both SQ and TR surfaces have a very small pinning

effect when W is small. The force needed to overcome the resistance from plane continues growing as the W increases. However, there is a remarkable dip down where W=[22,25] for SQ and W=[30,32.5] for TR surface, after that, the Boc continue growing, until the tile size is much larger

than the drop size, reaching a platform of maximum. 0.9 0.8

SQ TR

0.7 0.6

Boc

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

0.5 0.4 0.3 0.2 0.1 0.0 0

5

10

15

20

25

30

35

40

W Fig. 7 The critical Bo number for different tile size. To further gain insight into the mechanism that pins the drop, we inspected the surface force imposed on the drop. To this end, the coordinate is tilted as the plane for illustration clarity, and each pair force with a wall particle and a sessile drop liquid particle is summed up as the surface force. The surface force comprises three parts, the supporting force (z direction), the resistance force (x direction), and lateral force (y direction). The supporting force is verified having a mean value of the component of the gravity perpendicular to the plane. Only the downward resistance force is of great concern, and were plotted against the drop’s center position in Fig. 8. The value

ACS Paragon Plus Environment

19

Langmuir

is negative because the resistance force is against the drop advancing direction. Four cases with different pattern (SQ and TR) and different tile size (W=5 and W=10) are investigated, under the condition of 𝛼𝛼 = 25°. The drop starts to slide from X=0 point with zero velocity, after the initial irregularity, the velocity exhibits a periodic fluctuation with its period two times of W. In cases

with W=15 surfaces, we found the force curve has only one peak and one bottom in a cycle for SQ pattern, while two peaks and two bottoms for TR pattern. The surface force peaks at the similar position for both SQ and TR surfaces, while the positon where drop on TR surface experience the minimum surface force deviates from that on SQ surface. Moreover, it is noticeable that the peaks and bottoms for TR surfaces are milder than that in SQ surfaces. On the other hand, in cases with W=5 cases, the force curves are relatively flat, with no significant peaks and bottoms. -350

Sliding -300

Pinning

-250

Wall Force

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 28

-200 -150 -100

SQ W=5 TR W=5

SQ W=15 TR W=15

-50 0 0

15

30

45

60

75

90

105

120

135

150

X

Fig. 8 Force evolution during drop sliding. Blue arrows indicate two scenarios with different driving force. The analysis of the force evolution allows the understanding of the observed Boc difference for a

drop on different patterned surfaces. The Boc actually represents a static energetic barrier that must

be overcome by gravity before the drop start to move, in Fig. 8 the peak of the force curve

corresponds to the energetic barrier. To illustrate, the body force (driving force) is represented by

ACS Paragon Plus Environment

20

Page 21 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

blue arrows in Fig. 8. If the body force is lower than the peak, the drop simply pins (the lower arrow). Only if the body force exceeds the height of the peak, the drop will have continuous motion downward the plane (the upper arrow). Interpreting the serious of surface force, we could tell the Boc of SQ W=15 surface should be higher than that of TR W=15, the conclusion is supported by the Boc analysis in Fig. 7. The surface with W=5 have very flat surface force curve, indicating the Boc is even smaller, which is also proved by Fig. 7.

Although the surface force peak determines the minimum tilt angle to unpin the drop, the area under the force curve corresponds to the work of resistance force down along the drop sliding down. Comparing the four cases in Fig. 8, despite the peak and bottom value are quite different, the area under each curve are similar. This implies over a long term, the drop endures a similar averaging resistance which weakly depends on the pattern size and shape. This is insistent with what we found in Fig. 6 that the average sliding velocity is insensitive to pattern size and shape. 3.4

Wriggling movement

The stick slip phenomenon was reported on all chemical heterogeneous surfaces with sharp wettability interface 6. In this section we want to find out whether the center of drop move wriggling like the frontier. To this end, we make drops slide on different surfaces. The evolution of center’s velocity against the normalized position is plotted in Fig. 9. The normalized position is a relative position by W. A length of 4 times of W is displayed in Fig. 9. It is evident that, with W=5, the drops move smoothly on both SQ and TR surface, while drops on TR surface are faster. On W=10 surfaces, drops start to exhibit a velocity fluctuation with period of 2 times of W. The fluctuation is stronger in W=15 cases. On account of the energy dissipation of the obvious

ACS Paragon Plus Environment

21

Langmuir

fluctuation, the average velocity is lower than that on W=5. This is consistent with the finding that for a known Bo, the drop moves slightly slower if the W is bigger (see Fig. 6). 0.14 0.13 0.12 0.11

Velocity

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 28

0.10 0.09 0.08

SQ W=5 SQ W=10 SQ W=15 TR W=5 TR W=10 TR W=15

0.07 0.06 0.05 0.04 0

1

2

3

4

Relative position (W)

Fig. 9 Velocity of drop’s center on different surfaces. 3.5

Varied route patterns

Due to chemical pattern, drops may experience non-negligible unbalanced lateral force, which does not exist with homogeneous surface. In most cases, this lateral force only spurs the drop to swing slightly, however, in the following cases, the force is strong enough to lead to different advancing direction of the drop, which produces varied routes. Fig. 10 shows three main routes observed on TR W=20 surfaces with different tilting angles. We name the three route patterns, zigzag, oblique, and straight, respectively. In zigzag pattern, the drop make turns once it move forward a distance of W, leaving a zigzag trace. In oblique pattern, the drop is heading directly towards the tilted 30° directions with minor direction fluctuation,

ACS Paragon Plus Environment

22

Page 23 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

generating an oblique smooth aisle. In straight pattern, which is the most common pattern, the drop slides directly downwards with periodical shape deformation. The zigzag and straight routes are both bilateral symmetric, while the oblique route has two sub patterns: tilt to left or right. Fig. 10 only exhibits the right oblique trace, while the left is a mirror symmetry of that. By repeating the test for dozens of times, we found each side has even odds/each side has equal probability, and it is hard to predict the tilt direction from the beginning. Some initial disturbance is believed to be decisive for left or right. Once the oblique aisle is formed, further disturbance only leads to small fluctuation, not direction change. Besides the three main route patterns, more transition routes are found amid the zigzag/oblique and oblique/straight. When Bo falls between the typical value for zigzag and oblique, the route first advance in oblique direction, but make turns after random distance, making the trend of the drop quite unpredictable. We also observe shift between oblique and straight when Bo is higher than typical value for oblique but lower than the straight value. The varied routes are also found in SQ patterned surfaces. The details are not addressed herein. These remarkable phenomena caused by lateral force on chemical chessboard-patterned surface needs further investigation. It is thrilling if we are able to change the drop moving direction only by tilt the plane more or less. Deep understanding of this complex effect may lead to fully control of drop motion in open surface.

ACS Paragon Plus Environment

23

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

a. Zigzag (tilt 14°)

b. Oblique (tilt 18°)

Page 24 of 28

c. Straight (tilt 24°)

Fig. 10 Varied route patterns of drops with same surface. (View perpendicular to the plane)

4

Conclusion

In this study, we provide an adequate scaling relationship between MDPD model and an 80% w/w glycerol/water solution. Sliding of drops on chemical chessboard-patterned surfaces formed by square or triangular tiles with different tile sizes are systematically studied. Various morphologies of wetting area are found for the same volume drop resting on horizontal plane with different tile size. The area of wetting is fluctuating as the tile size changes, while the fraction of hydrophilicity shows a general increasing trend as the tile size grows. Like stripe-patterned surface, the pinning phenomenon is observed on tilted plane, the minimum body force needed to depin the drop has a strong non-linear relationship with both the tile size and shape. Instead, once the drops start to move, the slope ∆Ca/∆Bo depends weakly on the actual surface pattern size or shape. In the dynamic area, the drops slide smoothly if the tile size is much smaller than the drop, and will wriggle if comparable or bigger. Furthermore, three main route patterns (zigzag, oblique, and

ACS Paragon Plus Environment

24

Page 25 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

straight) are recognized, and more transit routes are also observed. Our study could be of help to design devices for passive drop manipulation by suitable tailoring of the chemical pattern, with potential applications in microfluidics and drop sorting.

5

Acknowledgements

This work was supported by the National Natural Science Foundation of China (Grant No. 51276130). The grant is gratefully acknowledged.

References (1)

Mistura, G.; Pierno, M. Drop Mobility on Chemically Heterogeneous and LubricantImpregnated Surfaces. Adv. Phys. X 2017, 2 (3), 591–607. (2) Morita, M.; Koga, T.; Otsuka, H.; Takahara, A. Macroscopic-Wetting Anisotropy on the Line-Patterned Surface of Fluoroalkylsilane Monolayers. Langmuir 2005, 21 (3), 911–918. (3) Jansen, H. P.; Bliznyuk, O.; Kooij, E. S.; Poelsema, B.; Zandvliet, H. J. W. Simulating Anisotropic Droplet Shapes on Chemically Striped Patterned Surfaces. Langmuir 2012, 28 (1), 499–505. (4) Suzuki, S.; Nakajima, A.; Tanaka, K.; Sakai, M.; Hashimoto, A.; Yoshida, N.; Kameshima, Y.; Okada, K. Sliding Behavior of Water Droplets on Line-Patterned Hydrophobic Surfaces. Appl. Surf. Sci. 2008, 254 (6), 1797–1805. (5) Sbragaglia, M.; Biferale, L.; Amati, G.; Varagnolo, S.; Ferraro, D.; Mistura, G.; Pierno, M. Sliding Drops across Alternating Hydrophobic and Hydrophilic Stripes. Phys. Rev. E 2014, 89 (1), 012406. (6) Varagnolo, S.; Schiocchet, V.; Ferraro, D.; Pierno, M.; Mistura, G.; Sbragaglia, M.; Gupta, A.; Amati, G. Tuning Drop Motion by Chemical Patterning of Surfaces. Langmuir 2014, 30 (9), 2401–2409. (7) Bonn, D.; Eggers, J.; Indekeu, J.; Meunier, J.; Rolley, E. Wetting and Spreading. Rev. Mod. Phys. 2009, 81 (2), 739–805. (8) Cristini, V.; Tan, Y.-C. Theory and Numerical Simulation of Droplet Dynamics in Complex Flows—a Review. Lab Chip 2004, 4 (4), 257–264. (9) Snoeijer, J. H.; Andreotti, B. Moving Contact Lines: Scales, Regimes, and Dynamical Transitions. Annu. Rev. Fluid Mech. 2013, 45 (1), 269–292. (10) Sui, Y.; Ding, H.; Spelt, P. D. M. Numerical Simulations of Flows with Moving Contact Lines. Annu. Rev. Fluid Mech. 2014, 46 (1), 97–119.

ACS Paragon Plus Environment

25

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 28

(11) Semprebon, C.; Varagnolo, S.; Filippi, D.; Perlini, L.; Pierno, M.; Brinkmann, M.; Mistura, G. Deviation of Sliding Drops at a Chemical Step. Soft Matter 2016, 12 (40), 8268–8273. (12) Ding, H.; Spelt, P. D. M. Wetting Condition in Diffuse Interface Simulations of Contact Line Motion. Phys. Rev. E 2007, 75 (4), 046708. (13) Borcia, R.; Borcia, I. D.; Bestehorn, M. Drops on an Arbitrarily Wetting Substrate: A Phase Field Description. Phys. Rev. E 2008, 78 (6), 066307. (14) Lohse, D.; Zhang, X. Surface Nanobubbles and Nanodroplets. Rev. Mod. Phys. 2015, 87 (3), 981–1035. (15) Li, Z.; Bian, X.; Li, X.; Deng, M.; Tang, Y.-H.; Caswell, B.; Karniadakis, G. E. Dissipative Particle Dynamics: Foundation, Evolution, Implementation, and Applications. In Particles in Flows; Advances in Mathematical Fluid Mechanics; Birkhäuser, Cham, 2017; pp 255– 326. (16) Lei, H.; Caswell, B.; Karniadakis, G. E. Direct Construction of Mesoscopic Models from Microscopic Simulations. Phys. Rev. E 2010, 81 (2), 026704. (17) Louis, A. A.; Bolhuis, P. G.; Hansen, J. P. Mean-Field Fluid Behavior of the Gaussian Core Model. Phys. Rev. E 2000, 62 (6), 7961–7972. (18) Warren, P. B. Vapor-Liquid Coexistence in Many-Body Dissipative Particle Dynamics. Phys. Rev. E 2003, 68 (6), 066702. (19) Chen, C.; Gao, C.; Zhuang, L.; Li, X.; Wu, P.; Dong, J.; Lu, J. A Many-Body Dissipative Particle Dynamics Study of Spontaneous Capillary Imbibition and Drainage. Langmuir 2010, 26 (12), 9533–9538. (20) Li, Z.; Hu, G.-H.; Wang, Z.-L.; Ma, Y.-B.; Zhou, Z.-W. Three Dimensional Flow Structures in a Moving Droplet on Substrate: A Dissipative Particle Dynamics Study. Phys. Fluids 2013, 25 (7), 072103. (21) Warren, P. B. Hydrodynamic Bubble Coarsening in Off-Critical Vapor-Liquid Phase Separation. Phys. Rev. Lett. 2001, 87 (22), 225702. (22) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117 (1), 1–19. (23) Anderson, R. L.; Bray, D. J.; Ferrante, A. S.; Noro, M. G.; Stott, I. P.; Warren, P. B. Dissipative Particle Dynamics: Systematic Parametrization Using Water-Octanol Partition Coefficients. J. Chem. Phys. 2017, 147 (9), 094503. (24) Pagonabarraga, I.; Frenkel, D. Dissipative Particle Dynamics for Interacting Systems. J. Chem. Phys. 2001, 115 (11), 5015–5026. (25) Groot, R. D.; Warren, P. B. Dissipative Particle Dynamics: Bridging the Gap between Atomistic and Mesoscopic Simulation. J. Chem. Phys. 1997, 107 (11), 4423–4435. (26) Li, Z.; Bian, X.; Tang, Y.-H.; Karniadakis, G. E. A Dissipative Particle Dynamics Method for Arbitrarily Complex Geometries. J. Comput. Phys. 2018, 355, 534–547. (27) Varagnolo, S.; Ferraro, D.; Fantinel, P.; Pierno, M.; Mistura, G.; Amati, G.; Biferale, L.; Sbragaglia, M. Stick-Slip Sliding of Water Drops on Chemically Heterogeneous Surfaces. Phys. Rev. Lett. 2013, 111 (6), 066101. (28) Arienti, M.; Pan, W.; Li, X.; Karniadakis, G. Many-Body Dissipative Particle Dynamics Simulation of Liquid/Vapor and Liquid/Solid Interactions. J. Chem. Phys. 2011, 134 (20), 204114. (29) Wang, Y.; Chen, S. Numerical Study on Droplet Sliding across Micropillars. Langmuir 2015, 31 (16), 4673–4677. (30) Wang, Y.; Chen, S.; Liu, Y. Spontaneous Uptake of Droplets into Non-Wetting Capillaries.

ACS Paragon Plus Environment

26

Page 27 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Langmuir

(31) (32) (33)

(34)

Comput. Fluids 2016, 134–135, 190–195. Chen, C.; Zhuang, L.; Li, X.; Dong, J.; Lu, J. A Many-Body Dissipative Particle Dynamics Study of Forced Water–Oil Displacement in Capillary. Langmuir 2012, 28 (2), 1330–1336. Kusumaatmaja, H.; Yeomans, J. M. Modeling Contact Angle Hysteresis on Chemically Patterned and Superhydrophobic Surfaces. Langmuir 2007, 23 (11), 6019–6032. Jansen, H. P.; Sotthewes, K.; Ganser, C.; Zandvliet, H. J. W.; Teichert, C.; Kooij, E. S. Shape of Picoliter Droplets on Chemically Striped Patterned Substrates. Langmuir 2014, 30 (39), 11574–11581. Kusumaatmaja, H.; Leopoldes, J.; Dupuis, A.; Yeomans, J. M. Drop Dynamics on Chemically Patterned Surfaces. Europhys. Lett. EPL 2006, 73 (5), 740–746.

ACS Paragon Plus Environment

27

Langmuir 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 28

Table of Contents Graphic

ACS Paragon Plus Environment

28