{'fraction_non_alphanumeric': 0.029315512074132764, 'fraction_numerical': 0.011092355919942128, 'mean_word_length': 4.6281504459092675, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Hydrodynamics of the Developing Region in Hydrophobic Microchannels: A Dissipative Particle Dynamics Study

Dissipative Particle Dynamics (DPD) is becoming a popular particle based method to study flow through microchannels due to the ease with which the presence of biological cells or DNA chains can be modeled. Many Lab-On-Chip (LOC) devices require the ability to manipulate the transport of cells or DNA chains in the fluid flow. Microchannel surfaces coated with combinations of hydrophilic and hydrophobic materials have been found useful for this purpose. In this work, we have numerically studied the hydrodynamics of a steady nonuniform developing flow between two infinite parallel plates with hydrophilic and hydrophobic surfaces using DPD for the first time. The hydrophobic and hydrophilic surfaces were modeled using partial-slip and no-slip boundary conditions respectively in the simulations. We also propose a new method to model the inflow and outflow boundaries for the DPD simulations. The simulation results of the developing flow match analytical solutions from continuum theory for no-slip and partial-slip surfaces to good accord. The entrance region constitutes a considerable fraction of the channel length in miniaturized devices. Thus it is desirable for the length of the developing region to be short as most microfluidic devices such as cell or DNA separators and mixers are designed for the developed flow field. We studied the effect of a hydrophilic strip near the inlet of a microchannel on the effective developing length. We find that the presence of the hydrophobic strip significantly reduces the developing length.

transport of cells or DNA chains in the fluid flow. Microchannel surfaces coated with combinations of hydrophilic and hydrophobic materials have been found useful for this purpose. In this work, we have numerically studied the hydrodynamics of a steady nonuniform developing flow between two infinite parallel plates with hydrophilic and hydrophobic surfaces using DPD for the first time.
The hydrophobic and hydrophilic surfaces were modeled using partial-slip and no-slip boundary conditions respectively in the simulations. We also propose a new method to model the inflow and outflow boundaries for the DPD simulations. The simulation results of the developing flow match analytical solutions from continuum theory for no-slip and partial-slip surfaces to good accord.
The entrance region constitutes a considerable fraction of the channel length in miniaturized devices. Thus it is desirable for the length of the developing region to be short as most microfluidic devices such as cell or DNA separators and mixers are designed for the developed flow field. We studied the effect of a hydrophilic strip near the inlet of a microchannel on the effective developing length. We find that the presence of the hydrophobic strip significantly reduces the developing length.

I. INTRODUCTION
Miniaturization of fluidic devices from bench-top to palm-top size has been progressed considerably in recent years [1]. The advent of novel micro-and nano-manufacturing and fabrication techniques has equipped scientists and engineers the ability to manipulate the transport of micro/nano liters of fluid through the micro/nano channels in these devices [2]. The applications of such microfluidic devices are spread through fields like electronicchip cooling, chemical synthesis, targeted cell isolation, bio-particle separation processes, chromatography, micro-particle sorting, micro-reaction, micro-mixing, genomic/proteomic studies and others [1,[3][4][5][6][7]. Shrinking the size of table-top labs has led to the evolution of new kinds of 'on-chip' bio-assays such as lab-on-chip, blood-on-chip, cell-on-chip and neurons-on-chip [1,5,7,8]. The advantages of these microfluidic devices are manifold: they are portable, fast, affordable, accurate and energy efficient. The small sample and reagent volumes required and ease of use by non-experts also equip them to cater to 'point-of-care' needs. Fluid flow through channels of a few microns in size is a common feature of all these devices.
Hydrophobicity arises when the surface energy of the solid-liquid interface is high and is usually measured from the contact angle made by a sessile drop on a surface [10,11]. Hydrophobic coatings on the walls of the microchannels facilitate larger flow rates compared to hydrophilic counterparts for the same pressure drop as they offer less resistance to flow [12].
Hydrophobic surfaces can also amplify electro-kinetic pumping, aid passive chaotic mixing and also mitigate the possibility of choking or adhesion of suspended analytes [13][14][15]. On hydrophobic surfaces, the traditional no-slip boundary condition is not valid and, instead, the fluid is modeled with a finite velocity at the wall. Navier proposed a generalized boundary condition to model the velocity of fluid (u) tangential to the wall, by assuming it to be proportional to the surface shear stress [16], where the slip-length β is the distance from the surface to the point where the linearly extrapolated velocity profile vanishes. The slip-length can be used to characterize the type of flow in channels; if β = 0 the flow is stick-flow (i.e. no slip), if β = ∞ the flow is plug-flow (e.g. shear free boundary) and any value of β between these two would represent partial slip-flow. In practical terms, the slip length at hydrophobic surfaces varies from a few nanometers to a few microns [17]. A slip length of up to 185 µm was reported in an experimental study [18], which is comparable to the size of boundary layers in macroscopic regime. In the macroscopic regime the no-slip boundary condition (i.e. velocity of the fluid at solid surface is equal to the wall velocity) captures the physics of flow adequately [9,19]. However, for flows in microchannels, partial-slip boundary conditions have to be applied at hydrophobic surfaces. In the remainder of this paper, we will use hydrophillic to refer to surfaces on which slip length is insignificant and the no-slip boundary condition is appropriate and hydrophobic for surfaces with a finite partial slip.
Numerical simulation of flow through microchannels is important for understanding the underlying physics of these flows, as well as minimizing effort and expense of experiments, especially during the design and optimization of microfluidic devices. The range of numerical simulation methods spans continuum based computational fluid dynamics (CFD) to atomistic level molecular dynamics (MD). Modeling flow through microchannels based on the continuum assumption suffers from several drawbacks. First, for flows with high Knudsen number (Kn) the continuum approximation may begin to fail [20]. Second, many of the flows include the presence of mesoscale particles such as DNA or individual biological cells. Treating the interaction of such second phase particles with the fluid medium and including the Brownian effects due to random thermal fluctuations becomes computationally prohibitive in CFD calculations. Microscopic modeling of above problems using MD is also not a practical choice as it is much too detailed and computationally expensive. As alternatives, discrete computational schemes like stochastic rotation dynamics, Brownian dynamics, lattice-Boltzmann method, smoothed particle dynamics, DPD have been developed primarily for the spatio-temporal scales relevant to these situations [21]. Among these discrete methods, DPD is a popular mesoscale scheme which bridges the gap between the macroscopic CFD and the microscopic MD [22].
DPD was introduced to study the dynamics of complex fluids such as colloids, soft matter and polymers [23]. A cluster of atoms or molecules are considered to form a single particle in the DPD scheme. The positions of the particles are updated in a Lagrangian framework, and thus DPD can be viewed as a coarse-grained version of MD. In the initial version of the method, the DPD particles were treated as point masses and larger sized particles The incorporation of appropriate wall boundary conditions for DPD has proved to be a challenge [21,25,26]. Discrete methods show density fluctuations near the walls. While such density fluctuations are realistic at the molecular level, they are considered spurious at the mesoscale continuum level at which DPD purports to model fluids. An instantaneous wall boundary (IWB) model introduced by Ranjith et al. [27] proved useful for modeling no-slip and partial-slip with minimum near wall fluid property perturbations. Inflow and outflow boundary conditions have also recently been modeled in DPD [28].
The hydrodynamics of the developing region inside channels with no-slip surfaces has been studied experimentally [29], analytically [30] and numerically using CFD [31]. Even though there have been considerable efforts to study the fluid transport in hydrophobic microchannels [9,17] in the fully developed regime, not much attention had been paid to the entrance effects with a few exceptions [32,33]. Since the entrance length is proportional to Reynolds number Re, it can reasonably be ignored in low Re flows. While typical microflows are characterized by Re < 30 [34], in a few microfluidic applications like micro heat-exchangers, micro-mixers etc the Re reaches the order of a few hundreds [15,35,36]. Wall shear stress effects and velocity distributions vary significantly at the entrance and these may eventually affect the separation efficiency of the microfluidic processes [32]. Moreover, the entrance region in hydrophobic channels is much longer than hydrophilic channels [33]. Reduction of entrance length of hydrophobic microchannels is very important for the design of some types of LOC devices. In the present work, we study the hydrodynamics of developing flow between two parallel plates with hydrophilic and hydrophobic surfaces using FDPD method for the first time. We also study the effect of hydrophilic patches in the entrance of the channel and their effect in reducing the developing length in microchannels.
In section II, analytical solutions for flow between two infinite parallel plates with partial slip are summarized briefly following [30,33,37]. In section III, the governing equations of the FDPD method are summarized. The details for the implementation of slip-wall boundary conditions and inflow and outflow boundary conditions are also presented. The hydrophobic walls with a short hydrophilic strip at inlet in section V. We find that the hydrophilic strip at the entrance shortens the development length significantly.

II. ANALYTICAL SOLUTION FROM CONTINUUM THEORY
In this section, we briefly summarize the analytical solution of developing flow between parallel flat plates following Sparrow et al [30], Chakraborty and Anand [33], and Duan and Muzychka [37] closely. Consider pressure driven flow between two infinitely long parallel plates separated by H = 2h. The steady, incompressible flow is governed by the mass and momentum balance equations in two dimensions, where u and v are velocities in x and y direction, ρ is the density and ν is the kinematic viscosity of fluid. The linearized momentum equation [30] is where u is the average cross sectional velocity. The analytical solution of the governing equation was obtained following Sparrow et al. [30] by assuming a no-slip boundary condition. Later, Chakraborty and Anand [33], and Duan and Muzychka [37] assumed that the fluid has a finite velocity at the wall and is modeled by the Navier boundary condition (Eq. 1) [33,37]. The slip-length β was modeled as a function of properties of the gas layer adjacent to the wall in case of flow between hydrophobic surface by [33] and as a function of Kn in case of fluid flow through microchannels for 0.001 < Kn < 0.1 by [37]. The governing equations are nondimensionalized using the hydraulic diameter D h (which is 2H for parallel plates), half width h and average velocity u. Hence the dimensionless parameters are taken form of Navier boundary condition is given by U = β ′ (dU/dη). The dimensionless steady velocity profile is a function of both spatial coordinates ξ and η and is given by [33] The eigenvalues α i satisfy .

III. SIMULATION PROCEDURE AND BOUNDARY CONDITIONS
A. Finite-size dissipative particle dynamics In this section, the formulation of FDPD model presented by Pan et al. [24] is summarized. The domain of interest consists of N DPD particles of finite size with a number density ρ. In this model, a set of molecules constitute an FDPD particle having a mass m i (i = 1, . . . , N) and mass moment of inertia I i (i = 1, . . . , N). The degree of coarse-graining depends on the degree of spatio-temporal detail required. Each FDPD particle obeys Newton's laws of motion and the translational motion is governed by the linear momentum equation where v i and f i are, respectively, the velocity of and the force on the ith particle. As the particles are of finite size, the angular momentum equation is enforced by where ω i is the angular velocity and f ij is the effective force exerted on the ith particle by the neighboring jth particle, at a distance r ij = r i − r j . The tangential forces are assumed to impart torques on the particles in proportion to the particle radii R i and thus . The position, linear, and angular velocities of each particle is determined by the total force exerted by the surrounding particles within a certain finite cut-off radius r c . In this scheme, the contribution from four types of forces (central (C), translational (T ), rotational (R) and stochastic (S)) are considered The total force on the ith particle due to the surrounding particles is given bỹ The total force on the ith particle may include external forces (f E ) from gravitational, magnetic or electro-osmotic forces, if any, The central conservative repulsive force acting along the line connecting centers is taken to be where a ij is the repulsion parameter, r ij = |r ij | andê ij = r ij /r ij is a unit vector. An appropriate weight function Γ(r ij ) is selected such that the conservative force decreases monotonically to 0 at r ij = r c . Most DPD simulations employ the form of Γ(r ij ) given by The translational force is assumed to have central and non-central dissipative components given by where γ ij C and γ ij S are the central and shear dissipation coefficients respectively. This frictional force attempts to reduce the relative velocity v ij = v i − v j between particles in both directions. The rotational dissipative force is taken to be of the form Finally, a stochastic force is also accounted through the expression given below where ∆t is the time step, d = 2 for two-dimensional simulations and tr[dW ij ] is the trace of symmetric independent Wiener increment matrix dW ij while dW ij A is its antisymmetric part. According to the fluctuation dissipation theorem, the random and dissipation coefficients are related by σ ij C = 2k B T γ ij C and σ ij S = 2k B T γ ij S . The stochastic forces and the dissipation forces together act to maintain a constant temperature during the simulations.
As mentioned in the Introduction, DPD simulations have shown spurious density fluctuations at walls when enforcing no-slip boundary conditions [25,26]. Recently, a new method of enforcing wall boundary conditions in FDPD simulations has shown substantial reduction in the density fluctuations [27]. In this method, when a fluid particle is within the range of influence of the wall, the particle interacts with the closest point on the wall as if there were a wall particle for that time step. The interaction of the wall particle and the fluid particle is separately specified. This proved to be a simple method to reduce spurious density variations as well as control slip at the wall. This method, referred to as the instantaneous wall particle boundary (IWB), was shown to be a computationally efficient procedure for modeling impenetrable walls. The slip velocity at wall was tuned by controlling the lateral dissipative force component between fluid and wall along the direction tangential to the wall. This is achieved by changing the lateral dissipation coefficient γ S pw = α(1 − r pw /r c ) 2 γ pp S , which acts only within a distance r c from the wall. Here r pw is the distance between wall and fluid particle. As in the fluid particle-particle interactions, the dissipative and random coefficients are related by σ pw S = 2k B T γ pw S . The slip-length increases as the slip modification parameter α is decreased. Thus boundary conditions ranging from no-slip to a large partial slip could be achieved by tuning a slip modification factor α described in that scheme [27].
We note that, some surfaces achieve superhydrophobicity by enhanced roughness which trap air pockets. In our simulations, we do not model roughness or the second phase fluids at the wall. Instead we specify an effective slip velocity using the Navier condition given by Eq. (1).

C. Determination of slip-length of a hydrophobic surface
A simulation box of 20r c × 30r c size taken to be periodic in the stream-wise direction and bound by IWB walls in y direction is used to estimate the slip-length β. For a constant body force (f E = 0.01), the volume flow-rate per unit area Q is calculated. The β ′ is estimated from the theoretical expression [38], where slip flow-rate Q slip for different hydrophobic surfaces were obtained by tuning the parameter α. As seen from Fig. 1, the flow-rate increases with increasing β ′ . Thus the flow-rate obtained for a certain applied body force increases with increasing partial slip in accordance with the experimental findings reported earlier [17,38].
The expression for fully developed velocity profile for slip flow in dimensionless form is given by [33], is 2.25 and that from FDPD simulation is 2.251. Moreover for slip flow (β = 2.25r c ) the flow-rate ratio from the simulation Q 30rc slip /Q 20rc slip is 1.954 while theoretical prediction is 1.95. The FDPD scheme along with the IWB wall is thus able to capture the developed flow hydrodynamics of hydrophilic and hydrophobic parallel plates.

D. Channel inflow and outflow conditions
In order to compare the FDPD simulations with analytical results, the inflow has been ensured to be uniform at the inlet with a velocity ( u = 1) at unit temperature (k B T = 1). To  Fig. 3). The repulsive interaction forces exerted by these particles on the fluid particles in outlet is calculated by trial and error. The FDPD velocity profiles were found to match the analytical fully developed profiles to good accord for an inter-particle conservative force parameter value of a po = κa pp with κ = 1.2.

IV. FDPD SIMULATIONS OF DEVELOPING FLOW
We devote this section to study the ability of the FDPD model for simulating developing flows in channels with no-slip and partial slip boundary conditions. We consider a steady developing flow between two long parallel plates. The size of the 2D simulation domain is taken to be 20r c × 400r c and filled with ρ = 3 particles per unit volume. The effective fluid viscosity is calculated to be µ = 1 [27]. The short range interactions were calculated with a cut-off radius r c = 1. The Reynolds number Re = (ρuH/µ) is calculated to be 60 for the geometry and fluid properties under consideration. The wall and fluid particles are taken to be of the same size and thus the particle size coefficients are λ ij = λ ji = 1 2 . The maximum repulsion parameter between fluid-fluid a pp = 75k B T /ρ is chosen according to Groot and Warren [22] and fluid-wall a pw = 20 is taken ('p' and 'w' represent fluid and wall particles respectively) according to Ranjith et al. [27]. The central and shear coefficients of the dissipative and random forces are taken to be γ C pp = γ S pp = γ C pw = 4.5 and σ C pp = σ S pp = σ C pw = 3. The time step was taken to be dt = 0.01. All particles in the domain are arranged randomly with zero initial velocity. A body force f E is assumed to act on each particle, which accelerates through the domain. The net momentum (U) of domain increases from 0 to a uniform value of 1. There domain is decomposed into 400×200 bins in x and y direction across the length and breadth of channel to obtain statistical averages of the velocity inside the domain. The component of velocity in the direction of flow was averaged over 2 × 10 5 iterations to get statistically accurate results.
The force f E is adjusted for each partial-slip boundary condition to maintain the flow rate of Q = 1. For the range of forces applied in the present simulations, a value of κ = 1.2 (ratio of the interaction coefficient of the fluid particle and inflow and outflow boundary particles and inter-particle interaction coefficients) ensured that Q = 1 and the velocity profiles are very close to the analytical solution. It was observed that for κ < 1.15 the particles close to the outlet accelerate and for κ > 1.25 they decelerate. In both cases the simulated velocity field did not match the analytical solution.
A. Flow in a long hydrophilic channel A uniform velocity profile at inlet with average velocity u at the inlet transforms to a parabolic velocity profile at the outlet with a maximum velocity 1.5u. Within the developing region the velocity is a function of both x and y. When the flow is fully developed, velocity profile is given by Eq. (17) and remains same further downstream.
The no-slip condition at the wall is obtained by modeling the solid boundary with a slip modification factor α = 3 as reported in our earlier work [27]. The inflow and outflow boundary conditions were implemented as discussed in section III D. The velocity profiles at different axial positions (ξ = constant) are plotted in Fig. 4 (a), along with the analytical solution. In Fig. 4 (b), the velocity profiles at different heights (η = constant) are extracted and compared with the analytical solution. The FDPD simulation and analytical solutions were found to be in good agreement. The channel length was chosen to be L c = 20H to minimize end effects. Some end effects are apparent within H/4 from the outlet in the form of velocity fluctuations (less than 12% of the maximum velocity).

B. Flow in a long hydrophobic channel
The effective slip at the hydrophobic surface is obtained by choosing an appropriate parameter α as discussed in section III B. Due to the non-zero velocity at the wall, the acceleration of the central laminar core is less compared to that of the no-slip case. The velocity gradients produced by partial-slip wall are smaller than the no-slip walls. Due to this, the developing length L e for superhydrophobic microchannels is greater than in hydrophilic channels [33].
The velocity profiles obtained through the FDPD simulation is compared to the analytical solution given by Eq. (5). The theoretical and computational results are in good agreement as can be seen in Fig. 5. This simulation shows that the FDPD scheme, in combination with IWB wall, can capture the hydrodynamics of the steady non-uniform developing region of partial-slip flow accurately. The effect of a no-slip region at the entrance on the hydrodynamics of fluid flow between two hydrophobic surfaces is discussed in the next section.

V. ENTRANCE REGION WITH A HYDROPHILIC STRIP
The developing length or entrance length is the stream-wise distance from the inlet to the point at which the boundary layers formed on both the walls merge. In an engineering sense, this is quantified by the distance along the flow direction at which the centerline velocity reaches 99% of the maximum fully developed velocity. There are several empirical correlations for the entrance length as a function of the Reynolds number for no-slip channels [31,37]. For moderate Reynolds numbers, the developing length constitutes a considerable portion of a miniaturized LOC device, and hence a reduction in the entrance length of the microchannel is highly desirable. For partial-slip flow, the development of the flow field occurs over a much greater length compared to no-slip wall boundary conditions. This is because the acceleration of the central core for a no-slip wall is greater than that of a slipwall for the same Re. The greater shear stress at the wall for a hydrophilic surface enables the presence of the wall to be felt inside the domain over a shorter distance compared to their superhydrophobic counterparts. The latter have a smaller magnitude of the velocity gradient near the walls. We explore the hydrodynamics of flow between two long parallel superhydrophobic surfaces with a hydrophilic strip at the inlet in this section.
Experiments on the entrance hydrodynamics in microchannels with no-slip boundary conditions with aspect ratio H/W = 1 (W being the width of the channel) were recently reported by Ahmad and Hassan [29] for a hydraulic radius (D h = H) ranging from 100µm to 500µm over a range of Re numbers from 0.5 to 200. The entrance length of a hydrophilic microchannel with D h = 200µm and Re = 60 is interpolated from the empirical correlation obtained from their experimental data (Eq. (6) of Ref. [29]) is 943 µm. Using an analytical approach, Chakraborty and Anand [33] have presented a correlation that relates the developing length and Reynolds number of slip surfaces, Here, D h = H, for a given slip length (β ′ ). However, to the best of our knowledge there is no experimental data available till date on the hydrodynamics of entrance region of microchannels with hydrophobic surfaces. Hence, the above analytical solution in Eqn. 19  In this section the effect of introduction of a no-slip strip having length l s at the inlet, just before the hydrophobic surface (schematically shown in Fig. 6) is discussed. The velocity profiles obtained through the DPD simulations closely follows the analytical solution, so Eq. To study the influence of the no-slip strips, a number of DPD simulations were carried out by varying the initial hydrophilic strip length from 5r c to 25r c at Re = 60. As the stream-wise velocity gradient is maximum near the wall, the entrance length of a mixed hydrophilic-hydrophobic surface was found by monitoring the velocity near the wall. The transition between the hydrophilic and hydrophobic surface results in velocity fluctuations.
These fluctuations are minimum along the centerline (η = 0) and greatest at the wall.
So the developing length with a hydrophilic strip (L s e ) at the entrance, is estimated by determining the distance from the inlet to the point where the axial velocity for η = 0.9 becomes constant. The effect of the length of the inlet hydrophilic strip on the developing length is given in Table I. It was found that for a 10r c hydrophilic strip, the developing length reduces from 83r c to 30r c (see Fig. 7) and the percentage reduction in the developing length is Le−L s e Le ≈ 66%, although the portion of hydrophobic surface replaced by hydrophilic surface is only (l s /L e ) ≈ 1 8 ≈ 12.5%. Thus, the combination of hydrophilic-hydrophobic surfaces drastically reduced the developing length L e . The full development of velocity profile of such an arrangement takes place over a shorter distance than that of pure hydrophobic surfaces, as shown in Fig. 8. This finding is expected to be beneficial for the design of microfluidic devices and to optimize the size of LOC devices with hydrophobic channels. It is noteworthy that, from the analytical solution (Eq. 5) the central line velocity at the end of hydrophilic strip (x = 10r c ) is 1.244 and the fully developed centerline velocity of hydrophobic surface is 1.3. We infer that, the l s should be long enough to accelerate the central core of hydrophilic region nearer to the developed velocities of the hydrophobic surfaces. If the l s is too long the center region accelerates, and take a longer distance to decelerate and reach a developed state as shown in Fig. 8(b). Conversely, if the l s is short the acceleration of central core is not enough to reach the uniform fully developed state and may take a longer distance, as seen in Fig. 8(e).

VI. SUMMARY AND CONCLUSIONS
In this work, a modified DPD method has been shown to effectively capture the hydrodynamics of developing flows in microchannels with no-slip and partial-slip boundaries.
The simulations with no-slip and partial slip wall boundary conditions were shown to have excellent agreement with analytical results. A new method to model inflow boundary condition is proposed to obtain a uniform inlet velocity profile. Similarly, the outflow boundary conditions were modified to prevent the fluid from accelerating out of the domain. The presence of a small hydrophilic strip at the inlet of a hydrophobic microchannel was found to significantly reduce the development length of the remaining hydrophobic channel. This finding can potentially be used by the designers of LOC devices to optimize the size of the microfluidic devices.