Adhesion and Separation of Nanoparticles on Polymer-Grafted Porous

Article pubs.acs.org/Langmuir

Adhesion and Separation of Nanoparticles on Polymer-Grafted Porous Substrates Kolattukudy P. Santo,† Aleksey Vishnyakov,† Yefim Brun,‡ and Alexander V. Neimark*,† †

Department of Chemical and Biochemical Engineering, Rutgers University, Piscataway, New Jersey 08854, United States DuPont Central Research & Development, Wilmington, Delaware 19803, United States

‡

S Supporting Information *

ABSTRACT: This work explores interactions of functionalized nanoparticles (NP) with polymer brushes (PB) in a binary mixture of good and poor solvents. NP−PB systems are used in multiple applications, and we are particularly interested in the problem of chromatographic separation of NPs on polymer-grafted porous columns. This process involves NP flow through the pore channels with walls covered by PBs. NP−PB adhesion is governed by adsorption of polymer chains to NP surface and entropic repulsion caused by the polymer chain confinement between NP and the channel wall. Both factors depend on the solvent composition, variation of which causes contraction or expansion of PB. Using dissipative particle dynamics simulations in conjunction with the ghost tweezers free energy calculation technique, we examine the free energy landscapes of functionalized NPs within PB-grafted channels depending on the solvent composition at different PB grafting densities and polymer−solvent affinities. The free energy landscape determines the probability of NP location at a given distance to the surface, positions of equilibrium adhesion states, and the Henry constant that characterizes adsorption equilibrium and NP partitioning between the stationary phase of PB and mobile phase of flowing solvent. We analyze NP transport through a polymer-grafted channel and calculate the mean velocity and retention time of NP depending on the NP size and solvent composition. We find that, with the increase of the bad (poor) solvent fraction and respective PB contraction, NP separation exhibits a transition from the hydrodynamic size exclusion regime with larger NPs having shorter retention time to the adsorption regime with smaller NPs having shorter retention time. The observed reversal of the sequence of elution is reminiscent of the critical condition in polymer chromatography at which the retention time is molecular weight independent. This finding suggests the possibility of the existence of an analogous special regime in nanoparticle chromatography at which NPs with like surface properties elute together regardless of their size. The latter has important practical implications: NPs can be separated by surface chemistry rather than by their size employing the gradient mode of elution with controlled variation of solvent composition. (FFF),18,19 and size selective precipitation20 (SSP). These techniques are designed to separate NPs mainly according to their size and shape. However, many applications require NP separation by surface chemistry (that may be manifested, for example, by the degree of hydrophobicity or surface heterogeneity), especially when NPs are modified by specific functional groups (“ligands”).21−23 In polymer chromatography, separation by chemistry, e.g., chemical composition or microstructure of polymer chains, is achieved at the so-called critical point of adsorption (CPA).38 At CPA, the repulsive entropic effect of polymer chain confinement in the pores of solids substrate is compensated by the attractive enthalpic (adsorption) interaction.24 The balance of entropic and enthalpic factors results in molecular weight-independent elution of polymers. The corresponding separation regime,

1. INTRODUCTION Nanoparticles (NPs) are being used extensively in modern-day technologies in versatile roles, including various applications in diagnostics and therapeutics1,2 as carriers for drug delivery, in imaging as contrast agents,3 biosensing,4,5 electrodes of supercapacitors,6 light emitting diodes (LEDs),7,8 solar cells,9,10 and nanoelectronic devices. 11 Many of these applications involve interaction of nanoparticles with polymers or polymer brushes (PBs) grafted to solid substrates. PBs doped by NPs are used in the fabrication of nanocomposites, sensors and biomedical devices,12 colloidal stabilization, lubrication,13 and manipulation of nanoobjects on surfaces.14 Emerging applications of PB-grafted substrates in NP separation and purification15 attract special attention to NP− PB interactions. Traditionally, techniques for NP separation stem from those employed for characterization and separation of polymers, including size exclusion chromatography (SEC),16 hydrodynamic chromatography (HDC),17 field flow fractionation © 2017 American Chemical Society

Received: August 17, 2017 Revised: September 14, 2017 Published: September 15, 2017 1481

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

density, NP size, and solvent composition. We described two basic mechanisms of NP−PB interactions, NP adhesion at the PB exterior and NP immersion into PB with adhesion to substrate, and calculated the free energies of these adhesion states and the energy barriers between them. Somewhat similar observations were made by Hua et al.50 who modeled the distribution of NPs of different sizes in the vicinity of a semiflexible PB-grafted support. Under certain conditions, NP concentration indicated the existence of a potential well near the brush support. CGMD simulation of Nasrabad et al.51 considered NP influence on PB structure. The authors modeled PBs grafted onto planar surfaces and in cylindrical pores. Strong adsorption of polymer chains onto NPs caused PB collapse. At a low grafting density, a strongly nonuniform PB outer boundary was observed with polymer chains clustered around the NPs that were located close to the surface. This phenomenon is qualitatively similar to the behavior of low density PBs in bad solvent. PB became more uniform when the polymer grafting density and NP concentration increased. Cao et al.52 modeled NP diffusion within a cylindrical channel with PB-grafted walls using DPD framework with Lennard-Jones quasiparticles. The authors varied the grafting density and the attraction between NPs and tethered chains and observed a transition from static NPs trapped by the chains to fast transport of NPs in the channel center. In recent works from our group, we investigated the PB conformational transition and NP adhesion to PB in a binary solvent at varying solvent composition49,53 using DPD simulations. The original ghost tweezers (GT) method was employed to calculate the free energy landscape of NP−PB interactions and analyze the equilibrium adhesion states and the energy barriers separating these states. The DPD model and the GT method proposed in refs 49 and 53 constitute the methodological foundation of this work. The rest of the paper is structured as follows. Section 2 describes the methodological framework of our study: details of DPD simulations, parametrization of the coarse-grained model of polymer, solvent and nanoparticles, methodology of calculations of the free energy landscape of NPs near PBs by the ghost tweezers method, calculation of Henry constant of NP adhesion to PB in slitlike channels with PB-covered walls, and modeling of NP transport through PB-grafted slit-shaped channels. In section 3, we present the simulation results and discuss the dependence of NP−PB interactions on the NP size and solvent composition, variation of the Henry constant, and the retention time with the reduction in the solvent quality during NP flow through PB-grafted channels and demonstrate and quantify the observed transition from the size-exclusion to the adsorption regime of NP separation, which is accompanied by the reversal of the sequence of NP elution. Conclusions are summarized in section 4.

often called liquid chromatography at critical conditions (LCCC), has been applied to separate polymers of different architectures,24 such as diblock and triblock copolymers, star polymers, and combs at the isocratic and gradient elution modes.25−27 CPA can be achieved by varying the solvent composition, which controls the solvent quality.28 In a good solvent, attractive polymer−substrate interactions are suppressed and steric interactions dominate, favoring separation by hydrodynamic size (size exclusion chromatography) where large molecules elute faster than their smaller counterparts. A similar order of elution occurs in hydrodynamic chromatography in capillary channels or packed columns.17 In the other limiting case of strong adsorption achieved with a poorer solvent, elution occurs in the reverse order with larger molecules retained longer. CPA is achieved at an intermediate composition when retention becomes independent of the molecular weight. Recently, we demonstrated that the critical conditions of polymer adsorption exist in the chromatographic columns packed with nonporous particles in which the separation occurs in the interstitial volume between particles (flow through channels).29,30 The present work studies the interactions between NPs and PBs in a binary solvent by means of dissipative particle dynamics (DPD) simulations. We explore the dependence of the free energy landscape of NPs in the vicinity of PBs on the particle size and solvent composition and predict NP transport and retention in PB-grafted pore channels. Our motivation is to determine whether the reversal of the elution order and CPA, found in polymer chromatography upon varying the solvent composition, can also be achieved in NP chromatography, thus leading to size-independent NP elution. Although this particular task has (to the best of our knowledge) never been tackled in the literature, the NP−PB interactions in general have been actively explored both experimentally and theoretically. Of our particular interest are theoretical studies of the free energy of NP adhesion and sorption by PB and its influence on NP mobility. NP−polymer interactions are complex and involve multiple scales that require theoretical approaches31,32 and computer simulations33−44 on a coarse-grained supramolecular level. Theoretically, NP−polymer interactions were explored using self-consistent field theory (SCFT),31,33,35−37,39 Brownian dynamics (BD),40 and dissipative particle dynamics (DPD) 38 simulations. Milchev et al. 45 compared the interactions of NPs with PBs and polymeric solutions. They noted that free energy of NP rises steeply as NP penetrates deeper into an expanded PB, which could not be explained solely by the dependence of polymer density on the distance to the substrate.46,47 Halperin et al.31 explored the dependence of the disjoining pressure between NP and PB on solvent quality using SCFT. The reduction of the solvent quality, which leads to the eventual collapse of PB, reduces the associated free energy barrier and favors NP penetration and absorption within PB. Merlitz et al. in their MD study48 observed that chain polydispersity reduces the free energy of NP immersion into PB. Zhang et al.40 considered the effect of attraction between NP and polymer using BD simulations. For larger NPs, the free energy has a minimum with respect to NP location, thus indicating a preferential depth of NP immersion into PB. Our group49 explored NP interactions with PBs using DPD and focused at the interplay between attractive forces caused by polymer adsorption at the NP surface and NP attraction to the substrate and entropic repulsive forces at varying grafting

2. MODELS AND METHODS 2.1. Dissipative Particle Dynamics Simulations. In this work, we employed standard DPD implementation: molecules are dissected into soft quasiparticles (beads) each of which represents a group of atoms together. The beads interact with effective potentials parametrized to reproduce targeted structural and thermodynamic macroscopic properties of the material. The beads move according to Newton’s equations of motion while acted upon by a force, fi =

(B) (D) (R) ∑ (F(C) ij + Fij + Fij + Fij ) i≠j

1482

(1) DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir where f i is the total force on the ith particle resulting from interactions with its neighbors. f i sums up the short-range conservative repulsion (B) (D) (R) forces F(C) ij , bond forces Fij , drag forces Fij , and random forces Fij . (C) Fij is the force arising from a harmonic repulsive potential between

(

overlapping beads: F(C) ij = aij 1 −

rij Rc

)r̂ for r ≤ R and zero beyond ij

ij

c

the effective bead diameter RC. rij = |ri − rj|, and aij is known as the conservative repulsion parameter. Beads connected by bonds interact via harmonic bond forces F(B) ij = KB(rij − re)r̂ij, where KB is the effective bond rigidity and re is the equilibrium bond length. The drag and random forces are implemented to maintain the Langevin thermostat54,55 (D) (R) (R) F(D) ij = − γw (rij)(riĵ . vij)riĵ ; Fij = σw (rij)θij riĵ

(2)

where γ is the friction coefficient, σ2 = 2γkT, and vij = vi − vj is the relative velocity between ith and jth beads. The weight functions w(D)

(

and w(R) are related as w(D) = (w(R))2 with w(R)(rij) = 1 −

rij Rc

) for r

ij

≤ Rc and zero for rij > Rc. θij is a random variable with Gaussian statistics. In this work, we use the conventional implementation of DPD: all beads have the same effective diameter Rc and the same intracomponent repulsion parameter aII. The DPD model was parametrized in our previous works 49,53 and effectively mimics polyisoprene natural rubber (PINR) PB in benzene (good solvent)− acetone (bad solvent) binary solution. The bead size Rc is set to 0.71 nm. The bead density ρb equals 3R−3 c and aII = 42kBT/Rc to effectively reproduce the overall compressibility of the system under consideration.53,56 DPD simulations were performed in reduced units with Rc as the unit of length and kBT as the unit of energy using LAMMPS57 software. Configuration snapshots were created using the Visual Molecular Dynamics (VMD) program.58 2.2. System Setup. To calculate the free energy landscape of NP− PB interactions, we use the ghost tweezers (GT) method with a simulation setup similar to the one developed in our previous work53 (Figure 1). The coordinate system is chosen so that the grafted substrate surface is parallel to the XY plane and corresponds to z = 0. The PB is constructed by tethering linear polymer chains to the solid substrate at a given grafting density. Each chain consists of 90 beads of type P that corresponds to polyisoprene molecular weight of 6120 Da. The substrate is formed by an array of immobile S beads arranged in a cubic order at a very high density of 19Rc−3, which effectively makes the substrates strongly repulsive toward all polymer or solvent beads. We consider two grafting densities, Γ = 0.6 and Γ = 0.36 chains per nm2. PB is immersed in the binary mixture of good (G) and bad (B) solvent components. Solvent composition (expressed as volume fraction of good component G in the solvent xG = NG/(NG + NB), where N is the number of beads of a particular type) determines the PB conformation: the PB expands in the good solvent and contracts upon addition of the bad solvent. The parameters for polymer−solvent and solvent−solvent interactions are taken from refs 49 and 53 to effectively represent PINR in acetone (bad solvent, ΔaBP= aBP − aII = 7.0 kBT/Rc) and benzene (good solvent, ΔaGP = aGP − aII = 1.5kBT/ Rc). The solvents are miscible with each other, ΔaGB = 1.5 kBT/Rc. We also consider another system with the bad solvent component denoted as B2, which interacts with the polymer more favorably than B, ΔaB2P = 4.0kBT/Rc. NPs of spherical shape are composed of beads arranged on a hexagonal simple close-packed lattice and linked by strong harmonic bonds. NP is formed by two types of beads. The two outer layers of NPs consist of beads labeled N, which favorably interact with polymer. The inner layers of NP are formed by the “core” beads of type C that strongly repel all other beads in the system to ensure that NPs are impenetrable to solvent and polymer. Because of high bond rigidity, NP maintained an approximately spherical shape. To study the effect of NP size, we consider NPs of three different radii, 4Rc = 2.84 nm, 6Rc = 4.26 nm, and 8Rc = 5.68 nm, each consisting of M = 740, 2496, and 5917 beads, respectively. Short ligand chains are tethered to the NP

Figure 1. (a) Schematics of free energy calculations by the ghost tweezers method. Polymers are attached to a substrate to form a brush immersed in a binary solvent. NP is tied to its immobile “twin ghost tweezer” particle via harmonic bonds. The NP−GT displacement Δz gives the force exerted by PB on NP. (b) Representative snapshot of the PB−NP system used for free energy simulations. Colors: white beads, polymer; yellow, substrate; blue, ligands; and green, NP. The ghost particle beads are shown as transperent red beads. Solvent beads are not shown for clarity. (c) Representative snapshot of NP in a polymer-grafted channel under flow. The solvent velocity v(z) is approximated by a parabolic profile indicated by the arrows.

surface to mimic NP functionalization. Each ligand chain is composed of six beads of type L. Side chains are distributed uniformly over the NP surface. In all simulations reported in this work, the ligand grafting density is set to 0.39 nm−2, which corresponds to 40 ligands for RNP = 4Rc, 90 ligands for RNP = 6Rc, and 158 ligands for RNP = 8Rc. The polymer−polymer (P−P), NP−NP (N−N, N−C, and C−C), ligand−NP (L−N), and ligand−ligand (L−L) bonds have a force constant KB = 120kBT/R2c and an equilibrium bond length re = 0.8Rc (0.57 nm). All parameters are listed in Table 1. The boxes for simulations of NPs with diameters 8Rc, 12Rc and 16Rc have dimensions 40 × 40 × 50 R3c , 40 × 40 × 52 R3c , and 40 × 40 × 56 R3c , respectively. The extension in the Z dimension provides sufficient space for NP and PB to secure the uniformity of the bulk solvent outside PB. To analyze the dependence of the NP−PB adhesion on the NP size and solvent composition, we explore three characteristic systems differing by PB density and solvent−polymer interactions. In system 1, the PB consists of 484 polymer chains uniformly grafted on the substrate with the grafting density Γ = 0.6 nm−2. In system 2, the grafting density is reduced to Γ = 0.36 nm−2, which amounts to 289 polymers uniformly grafted on the substrate. It is worth noting that we considered a limited range of grafting densities (0.36−0.6 nm−2) at which PB is homogeneous. At lower grafting densities, one may expect additional effects related to the inhomogeneous distribution of grafted chains as the solvent quality worsens. These effects deserve a special investigation. In system 3, we improve the quality of the bad solvent to make bad solvent−polymer interactions more favorable; the grafting density is the same as in system 1, Γ = 0.6 nm−2. The bad solvent in system 3 is denoted as B2. In these systems, we obtain adhesion free energy landscapes for NPs of three sizes at five (systems 1 and 2) and six (system 3) different solvent compositions. For each system, we calculate the free energy landscape using the GT method explained below. Calculation of the free energy of NP at given distance z from the substrate requires an independent simulation of 1 million time steps or more. Altogether, we report the results of more than 590 independent simulations of at least 1 million time steps. 1483

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir Table 1. Parameters for Short-Range Repulsion and Harmonic Bond Potentialsa repulsion parameters aIJ, kBT/Rc I\J

C

N

L

A

P

G

B/B2

S

C N L A P G B/B2 S

42.0

60.0 42.0

60. 42.0 42.0

0.0 0.0 0.0 0.0

60.0 38.5 42.0 0.0 42.0

60.0 43.5 43.5 0.0 43.5 42.0

60.0 43.5 49.0/46.0 0.0 49.0/46.0 43.5 42.0

60.0 42.0 42.0 0.0 42.0 42.0 42.0 42.0

harmonic bond parameters type

KB (kBT/Rc2)

re/Rc

type

KB (kBT/Rc2)

re/Rc

P−P L−L

120 120

0.8 0.8

C−C GT

120 0.01

0.8 0

a

The repulsion parameter aij (in units kBT/Rc) between the beads and the bond parameters, the force constants KB, and the equilibrium bond length re used in the simulations. Bead types are denoted as C, core beads; N, NP surface beads; L, ligand beads; A, ghost particle beads; P, polymer beads; G, good solvent; B, bad solvent 1; B2, bad solvent 2; S, substrate. 2.3. Calculation of NP−PB Interaction Free Energy by the Ghost Tweezers Method. To measure the free energy of interactions between NP and PB, we apply the ghost tweezers (GT) method,53 which in silico mimics lab experiments with optical tweezers.59 The GT method is implemented by introducing a ghost particle that is an “identical twin” of the real NP (Figure 1a). The ghost twin (GT) particle is composed of beads of type A and of the same number of beads as the real NP kept in the undisturbed hexagonal order. GT is pinned at given position RGT; it is immobile and does not interact with any system components except for the NP. The NP is bound to GT by M harmonic springs of rigidity kGT connecting each respective pair of real and GT particle beads (i = 1, ..., M). Effectively, the NP is tethered to the given point RGT by an effective force exerted by the GT, FGT = ∑i M= 0kGT(Ri − RGT,i), where Ri is the fluctuating position of the ith bead of NP and RGT,i is the fixed position of its counterpart bead of the GT particle. This force counterbalances the sum of the forces between NP and polymer and solvent beads at the given NP position. Provided the system equilibration in the course of DPD simulation with the GT placed at RGT, we calculate the average distance z between the NP center of mass and substrate surface and the average GT force, FGT = ⟨FGT⟩, which due to the lateral uniformity of the system may have only one nonzero normal component that depends on the deviation, Δz = z − zGT, of the NP position z from the GT position zGT, FGT = KGT(z − zGT). Here, KGT = MkGT is the cumulative spring constant of the effective harmonic force acting between the real and ghost NPs. kGT is chosen to optimize the efficiency of calculations. At the beginning of a simulation, the GT is pinned to a certain location z0 far enough from PB to ensure a negligible NP−PB interaction so that at equilibrium, NP fluctuates around z0, Δz = 0 within the accuracy of the simulation and, respectively, the measured GT force, FGT = 0. Then, the GT is moved along the Z-axis in small increments toward the substrate. At each GT position zGT, the system is equilibrated allowing for the effective force acting between the NP and PB to be counterbalanced by the GT force, which is determined as a function of the NP position, FGT(z), from the average deviation of the NP center of mass from the tether point as described above. FGT(z) > 0 corresponds to an effective repulsion between PB and NP, whereas FGT(z) < 0 corresponds to effective attraction. At the equilibrium positions, FGT(z) = 0. With the simulation length we afforded in this work, an average error for the force data points is approximately 0.5kT/Rc, which is sufficient to determine the free energy landscapes that display different qualitative behaviors with reasonable accuracy. The change in the Helmholtz free energy, A(z), is determined as the work done by the GT in the course of quasi-equilibrium pulling of NP from the bulk solvent (position z0) to position z at the substrate,

A(z) =

∫z

z

FGT(z) dz

(3)

0

The NP location z0 far from the substrate serves as a reference point where A = 0. The free energy landscape, A(z), allows one to estimate the equilibrium distribution of NP near the substrate, determine the excess Gibbs adsorption of NPs, and calculate the Henry constant that quantifies the adhesion equilibrium. Let us consider a slitlike pore of half-width w with walls grafted with PB. The probability density P(z) of the NP location can be calculated as P(z) =

e−A(z)/ kBT 2w

∫0 e−A(z)/ kBT dz

(4)

The partition coefficient K of the nanoparticle between the polymer brush and the bulk solvent can be obtained using the Gibbs excess adsorption theory and by calculating the Henry constant, as detailed in the previous works.29,30,60 The Henry constant KH is defined from the excess adsorption of NPs per unit pore surface area, which represents the difference between the total amount of NPs in the channel Ntot and the amount of NPs in a reference volume of bulk fluid with the given concentration of NPs denoted as cB. This concentration is considered as low so that interactions between NPs are ignored. Because our ultimate goal is to understand the specifics of NP separation in pore channels with PB-grafted pore walls, it is reasonable to define the reference bulk volume as the volume available for the solvent flow. As shown below, the solvent flow does not penetrate into the PB, and the velocity profile is approximately parabolic (Figure 1c), which would be observed in a channel with solid walls of half-width wH = w − wPB. We define wH as the hydrodynamic width of the PB-grafted channel and wPB as the hydrodynamic thickness of the PB. The reference bulk volume for the definition of excess adsorption Nex represents the volume of the mobile phase outside of the PB

Nex =

1 Ntot − c BwH = c B 2

∫0

w

e−A(z)/ kBT dz − c BwH

(5)

Here, Nex is the total amount of NPs per unit substrate surface area. Accordingly, the Henry constant is given by KH =

Nex = cB

= wPB+

∫0

∫0 w

w

e−A(z)/ kBT dz − (w − wPB)

(e−A(z)/ kBT − 1) dz

(6)

As the reference technique, we consider the hydrodynamic chromatography employed for size exclusion separation of polymers 1484

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir and NPs alike.17,61,62 In HDC, the separation is achieved by the particle size due to the inhomogeneity of the solvent flow in the pore channels with solid walls. It is worth noting that in this case the Henry constant defined by eq 6 in the absence of PB (wPB = 0) is negative and equal to the NP radius RNP, KH = −RNP, which reflects NP depletion (negative excess adsorption) at the solid surface (A(z) = ∞ at z < RNP), which is characteristic for HDC. 2.4. Modeling NP Transport. To model NP transport through a chromatographic channel, we simulate the solvent flow through a slitlike channel with PB-grafted walls. The system is periodic in the lateral directions. We consider 16 systems, which differ by the solvent composition and PB grafting density corresponding to systems 1−3. The simulation boxes in solvent flow modeling have the size of 40 × 40 × 70 R3c (28.4 × 28.4 × 49.7 nm3). Selected simulations of NP flow are performed in a larger simulation box of size 40 × 40 × 90 R3c (28.4 × 28.4 × 63.9 nm3). We first equilibrate the systems at a given solvent composition and then apply a constant force in the X direction on the solvent particles and run the simulation until steady flow is reached. The velocity profile vs(z) of the solvent across the channel is obtained by averaging the solvent particle velocities over a long time (100000 steps). For all systems considered, the velocity of solvent particles inside PB is zero, whereas in the channel outside PB, the velocity profile tends to be parabolic similar to the Poiseuille flow in channels with solid walls, see Figure 1c. This observation allows us to introduce the hydrodynamic thickness wPB of PB as the cutoff of the parabolic approximation assuming that solvent velocity is zero at z < wPB and

(

which follows the Poiseuille law, vsPA(z) = vm 1 −

(w − z)2 (w − wPB)2

and zero for other values of z. The solvent density is constant outside the brush, ρs(z) = ρb (ρb = 3R−3 c ), whereas inside the brush, ρs(z) = ρb − ρPB(z), where ρPB(z) is the PB density. In estimating the mean NP velocity, we followed the assumption61,63,64 accepted in the conventional HDC with solid wall channels that the velocity of particles positioned at distance z from the wall equals the solvent velocity at this distance. The mean velocity of NPs is determined by the distribution of solvent velocity across the channel. (z) is assumed to be Poiseuillian and The fluid velocity in HDC, vHDC n is given by eq 9 with wPB = 0 and χ = 2/3. Provided that the particle center cannot be closer to the wall than its radius RNP, the mean particle velocity reduced to the mean solvent velocity equals

τpHDC τs

), at w

PB

⟨vp⟩ =

⟨vsPA ⟩ =

(

τp τs

w

(w − z)2

∫0 ρs (z) dz

) dz = v χ m

(8)

⟨vPA s ⟩/vm,

the mean solvent velocity normalized by the maximum χ= solvent velocity, defined by eq 8. The normalized velocity profile in the parabolic approximation is then

⎛ (w − z)2 ⎞ ⎟ at wPB ≤ z ≤ w vnPA(z) = χ −1 ⎜1 − (w − wPB)2 ⎠ ⎝

(10)

1 1

1 + δ − 3 δ2

≈1−δ (11)

∫w

PB + RNP

e−A(z)/ kBT vn(z) dz

w

∫0 e−A(z)/ kBT dz

(12)

=

1 ⟨vp⟩

(13)

In the case of poor adsorption of NPs to PB, one may expect that eq 13 converts to the HDC eq 11 with δ = RNP/(w − wPB). The parabolic approximation (eq 9) with the hydrodynamic thickness wPB of PB determined from the solvent flow simulation in one channel of a particular width allows one to estimate the mean NP velocity and, respectively, the retention time for the channels of any width. This is particularly useful considering the fact that simulations of NP flow through channels suffer from serious computational difficulties. First, the real flow velocities in chromatographic columns are very low compared to the thermal velocities of particles in DPD simulations. Thus, the force applied to create such a flow is extremely weak, and it is hardly possible to obtain necessary statistics on NP diffusion across the channel cross-section even in a reasonably long simulation. Second, simulations of large channels are computationally unfeasible. In fact, we have shown with selected simulations of NP flow that the assumption that the NP velocity can be approximated by the solvent velocity in the absence of NP does not cause any significant error (see Supporting Information (SI), section 1). It is worth noting that the parabolic approximation for NP velocity is a simplified model; it does not consider the effects of possible partial NP penetration into PB, which would introduce additional frictional forces to NP flow. It also ignores the NP rotational motion that is pronounced especially at the edge of PB.

(7)

(w − wPB)2

NP

1 2 δ 3

The ratio of the particle retention time τp to the solvent retention time τs is calculated as

Eq 7 accounts for an inhomogeneous distribution of the solvent in the channel due to its depletion within PB and the fact that solvent located within PB is effectively immobile at z < wPB. Note that the mean retention time of the solvent is the reciprocal of ⟨vs⟩. The normalized velocity profile vn(z) = vs(z)/⟨vs⟩ is found to be independent of the applied force. In the parabolic approximation, the mean velocity is given by w

=

w

w

0

vnHDC(z) dz = 1 + δ −

The latter approximate equality holds at δ ≪ 1. In the case of NP flow in a PB-grafted channel, we assume that NP velocity is equal to the solvent velocity given by eq 9 at wPB + RNP ≤ z ≤ w and equals 0 at z < wPB + RNP, i.e., when NP is partially or completely immersed into PB, it is assumed to be immobile. Because the probability of NP location at distance z is proportional to the Boltzmann factor e−A(z)/kBT, the mean NP velocity (reduced to the mean solvent velocity) is calculated as

w

vm ∫ ρs (z) 1 −

w

represents the ratio of the where δ = RNP/w. The reciprocal of to the solvent retention time τs particle retention time τHDC p

∫0 ρs (z)vs(z) dz ∫0 ρs (z) dz

∫R

⟨vHDC ⟩ p

≤ z ≤ w. vm is the maximum velocity of the flow in the channel center chosen from the condition of equality of the total solvent flux to that given by the parabolic approximation. This means that the solvent flow in a PB-grafted channel is similar to the Poiseuille flow in a solid wall channel of half-width w − wPB. Thus, defined hydrodynamic PB thickness depends on the solvent composition and decreases as the good solvent fraction xG decreases, causing PB contraction. The hydrodynamic thickness of the brush wPB is found to be independent of the applied driving force and, respectively, of the maximum velocity vm. As such, the hydrodynamic thickness of PB does not depend (within the accuracy of our simulations) on the channel width, which allows us to predict the flow velocity profile in the channels of different widths based on the wPB value determined in the simulation box of the chosen size. To compare different systems, the solvent flux across the channel, the mean solvent velocity, has to be kept constant. To secure this condition, the solvent velocity vs(z) is normalized by the mean solvent velocity,

⟨vs⟩ =

1 w − RNP

⟨vpHDC⟩ =

3. RESULTS AND DISCUSSION 3.1. Free Energy of NP Adhesion to the Polymer Brush. We study the dependence of the adhesion free energy on the solvent quality by varying the good solvent fraction xG from 0.99 to lower values up to 0.6. At xG = 0.99, PB is fully

(9) 1485

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

effective NP−PB force and the free energy landscape as functions of the NP coordinate (distance of NP center from substrate surface) for system 1 (G-B solvent mixture, Γ = 0.6 nm−2). At xG = 0.99, we observe no adhesion of NP regardless of its size: the effective force is zero in the channel center beyond PB and becomes repulsive as NP gets closer to PB (Figure 4a). Repulsion increases steeply as NP is immersed into PB as it was observed in previously published works.50,53 As solvent quality is reduced (xG decreases), and an enthalpic attraction between NP and PB emerges because interaction of NP surface (or ligands) with polymer is more favorable than with the bad solvent component. In the region of attraction, the force is negative; GT prevents NP from penetrating into PB. The position of zero force corresponds to the adhesion equilibrium: at smaller distances, entropic repulsion overpowers the attractive forces, which prevents the NP from further immersion into the PB. Already at xG = 0.92, system 1 exhibits slightly attractive force (negative values) at the PB surface, but repulsion increases steeply as the NP is submerged into PB. Attraction becomes progressively more pronounced as the bad solvent fraction increases (xG decreases) and PB contracts. The location of the minimum of the force is also shifted closer to substrate, which results not only from the enthalpic attraction itself but also from the decrease in the brush thickness (Figure 2). The lower panels in Figure 3 display the dependence of the Helmholtz free energy calculated by integrating the measured force using eq 3 on the position of the NP center of mass. For integration, simulated force data points are smoothened using moving averages (solid lines of the force profiles). At xG = 0.99, free energy landscapes are monotonic and do not indicate any particular equilibrium state (because NP−PB interactions are repulsive). This behavior is reminiscent of that in HDC. As the bad solvent fraction increases, free energy landscapes exhibit well-defined minima that correspond to the adhesion

expanded; as solvent quality worsens, PB gradually transforms to a collapsed state. This is shown in Figure 2, which displays

Figure 2. PB density (monomers per unit volume) profiles at different good solvent fractions in the simulations. Red symbols corresponds to system 1 simulations with polymer grafting density 0.6 nm−2 at xG = 0.99 (circles), 0.90 (asterisks), and 0.80 (squares). Green symbols represents densities in system 2 with grafting density 0.36 nm−2, where circles correspond to xG = 0.99, asterisks to 0.89, and squares to 0.80. Blue symbols corresponds to system 3 with bad solvent 2 and grafting density 0.6 nm−2. Here, densities at xG = 0.99 (circles), 0.72 (asterisks), and 0.60 (squares) are shown.

density profiles obtained in different simulation systems: as xG decreases, the brush becomes denser and its effective width decreases. For example, in system 1, PB at xG = 0.8 is effectively ∼15% thinner and, respectively, denser than PB at xG = 0.99. Solvent composition strongly affects the interaction between the brush and NP. Figure 3 displays the dependence of the

Figure 3. Force (upper panel) and free energy (lower panel) profiles (right y-axis) in system 1 at different good solvent fractions xG ranging from 0.80 to 0.99 for different NP diameters 8Rc (blue), 12Rc (green), and 16Rc (red). z-coordinate corresponds to the position of the NP center of mass. Points show the actual values of the effective force obtained with the GT method, and the curves display the smoothed data (moving average). Free energy curves are obtained by integrating the smoothed data. Note the different scales of the right y-axis chosen to show the magnitude of the respective values of force and energy. Regions of negative force correspond to the NP−PB attraction, and zero force values indicate free energy minima corresponding to the equilibrium adhesion positions. Polymer grafting density is 0.6 nm−2. Black squares show the local density profiles of PB beads (left y-axis). 1486

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

Figure 4. Characteristic equilibrium positions of NPs in size exclusion (upper panels) and adsorption (lower panels) regimes: Systems (a) 1, (b) 2, and (c) 3. The NP diameter is 16Rc. Colors: yellow, substrate; pink, polymer; green, NP; and blue, ligands.

Figure 5. Force (upper panel) and free energy (lower panel) profiles of the NPs in system 2 at different good solvent fractions. z is distance from the substrate surface expresed in Rc. The polymer grafting density of the brush is 0.36 nm−2. Colors: black, brush density; blue, green, and red are data with NPs with diameters of 8Rc, 12Rc, and 16Rc, respectively. The representation scheme is the same as in Figure 3.

adsorption interactions prevail due to increased NP penetration into PB because of their smaller size. As xG decreases (xG = 0.87 and 0.80), the adsorption energy increases with the NP size. The intermediate solvent compositions of xG= 0.9−0.92 correspond to a region of transition from the entropydominated to the enthalpy-dominated NP−PB interaction, or from the size exclusion mode to the adsorption mode of NP separation. As discussed in the Introduction, a similar behavior is typical for adsorption of polymers, where there exists a compensation point known as the critical point of adsorption (CPA), where the enthalpic attraction is counterbalanced by the entropic repulsion regardless the polymer size. For polymers, CPA is characterized by the size independence of the Henry constant and retention time. We could not identify such a point for NP adhesion; however, the observed transition from repulsive to adsorption regimes occurs in a narrow range of solvent compositions for all of the systems studied.

equilibrium with NP partially immersed into PB. The corresponding snapshots are shown in Figure 4a. The location of equilibrium and the corresponding equilibrium adhesion energy depends on the NP size in a nontrivial fashion. When the good solvent fraction is high, that is, when PB is fully expanded, repulsion is stronger for larger NP: the larger the NP, the more pronounced is the restriction on the conformational freedom of the polymer chains. When the good solvent fraction is low and PB is denser, the adsorption interactions prevail and larger NPs have larger adsorption energy due to larger NP surface area on which polymer segments can be adsorbed. Thus, there is a crossover upon the solvent quality variation in the adsorption free energy dependence on the NP size, which is demonstrated by the free energy profiles (Figure 3). At sufficiently high good solvent fraction, attraction is more pronounced for smaller NPs: for example, at xG = 0.92 and 0.90, smaller particles have wider regions of attraction where 1487

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

Figure 6. Force (upper panel) and free energy (lower panel) profiles of NPs in system 3. Polymer grafting density is 0.6 nm−2, which is the same as in system 1, but the bad solvent component is B2, which interacts more favorably with the polymer than the one in system 1. Color and representation schemes are the same as those in Figure 3.

Figure 7. (a) Snapshot of NP (RNP = 6Rc) at the closest proximity (z = WPB + RNP) to PB in the size-exclusion regime at xG = 0.99. The bad solvent fraction is negligible, and NP is repelled from the brush. (b) Snapshot of adsorbed NP on PB of density 0.6 nm−2 at xG = 0.80. Colors: yellow, bad solvent beads; green, NP; blue, ligands; transparent pink beads, polymer beads. Bad solvent is excluded from PB, causing its contraction and enhancing NP−PB interactions. (c) Solvent (red, bad solvent; green, good solvent) and PB (black) density profiles at xG = 0.80 (squares) and xG = 0.99 (circles).

adsorption modes are given in Figure 4c. A better interaction between bad solvent and polymer radically weakens the enthalpic attraction between NP and PB because the latter results from the contrast between NP−PB and solvent−PB interactions. In system 3, we observe the same crossover between the size exclusion mode, where repulsion is stronger for larger NPs, and the adsorption mode with the opposite trend of stronger attraction for larger NPs. However, because of the generally weaker attractive forces in system 3, the crossover region shifts to lower xG values compared to systems 1 and 2; in system 3, it occurs at xG ≈ 0.72 (Figure 6). Overall, we conclude that the transition between the size exclusion and adsorption modes is observed for all systems studied, and the conditions at which the crossover occurs depend on the PB density and specifics of solvent−polymer interactions. Adhesion of NP to PB results from the preferential interaction of the polymer with NP over the bad solvent. Bad solvent is excluded not only from the exterior of PB but also from the gap between PB and NP that causes deeper penetration of NP into PB and increased effective NP−PB attraction. The strength of NP−PB attraction therefore

Figure 5 shows the effective forces and free energy landscapes for system 2, which is similar to system 1 except for a lower grafting density (Γ = 0.36 nm−2). The PB in system 2 is effectively ∼4Rc (∼21%) thinner than in system 1 (Figure 2). Reduced PB density leads to less efficient entropic repulsion, and therefore, the interval of z values corresponding to attractive NP−PB interaction widens. This causes significantly stronger adhesion: free energy minima at lower xG in system 2 are much lower than those in system 1 at the same compositions. However, the transition from size-exclusion to adsorption modes occurs in the same range of compositions as in system 1 at xG = 0.89−0.92. Essentially, the crossover behavior is similar to that of system 1. Figure 4b shows snapshots of the equilibrium positions of the NP of diameter 6Rc with different solvent compositions. In system 3, the polymer grafting density is the same as in system 1 (Γ = 0.6 nm−2), but the bad component interacts with polymers more favorably (aP,B2 = 46.0kBT/Rc instead of 49.0kBT/Rc in systems 1 and 2). The resulting force and free energy profiles in system 3 are shown in Figure 6. The snapshots of NP equilibrium states in size exclusion and 1488

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

Figure 8. Solvent velocity and density profiles across PB-grafted channel of width 66.75Rc at different good solvent fractions xG in systems (a) 1, (b) 2, and (c) 3. Computed solvent velocity vs(z) (dashed lines), normalized solvent velocity vn(z) (solid lines), PB density (circles), and binary solvent density (crosses) are shown. Red, green, blue, cyan, and magenta colors represent xG= 0.80, 0.87, 0.90, 0.92, and 0.99, respectively, in (a and b) and 0.60, 0.70, 0.72, 0.80, and 0.99, respectively, in (c).

composition because the latter affects the brush extension and, respectively, the hydrodynamic width of the channel. Figure 8a presents the solvent velocity profiles at good solvent fractions 0.99, 0.92, 0.90, 0.87, and 0.80 in system 1 with PB density Γ = 0.6 nm−2. The channel width, which is the distance between the substrate surfaces in the Z direction is 66.75Rc (47.4 nm). The profiles show that the solvent velocity is zero inside PB, but within the space between the brushes, the solvent velocity profiles tend to be parabolic as in a Poiseuille’s flow between solid walls. The density distributions of the solvent and the polymer are correlated with the velocities in the outer region of the brush, where the velocity profile changes from zero to the parabolic type. As shown in Figure 8, the solvent density ρs(z) = ρb = 3Rc−3 in the bulk region (outside the brush), whereas inside the brush, ρs(z)= ρb − ρPB(z). Thus, there is a considerable amount of solvent inside PB, and this solvent is practically immobile. In HDC channels without PB, eq 8 integrates to 2/3vm (see section 2.3), and thus, the normalized velocities vn(z) have a peak of 3/2. The presence of PB decreases the mean velocity because a part of the solvent that is inside PB has zero velocity. Integrating eq 8, taking into account that the solvent density is constant in the bulk (z > wPB) and is approximately given by ρb − ⟨ρPB⟩ in the interior of the brush, we obtain

depends on the strength of repulsion between polymer and bad solvent. Panels a and b in Figure 7 show snapshots of NP−PB configurations at xG = 0.99 and 0.80 for PB of density 0.6 nm−2. Visual comparison of these two snapshots shows the specifics of NP−PB adhesion in binary solvents. In the absence of bad solvent (Figure 7a), PB prefers to be surrounded by good solvent, effectively repelling NP. Addition of bad solvent (Figure 7b), causes PB contraction and solvent exclusion between NP and PB that enhances NP penetration into PB. The differences in the polymer and solvent density profiles shown in Figure 7c quantify the change in the composition of the PB interface upon the addition of bad solvent. 3.2. Solvent Flow in the Polymer-Grafted Channel. Chromatographic separation of NPs depends not only on the adhesion interactions but also on the specifics of the solvent flow through polymer-grafted channels that are also controlled by the solvent quality. Figure 8 presents the results of simulation of solvent flow through slitlike channels between PB-grafted walls as shown in Figure 1c. We consider the same PB−binary solvent systems 1−3 as in the free energy calculations. The simulation boxes are open and periodic in the XY plane. A constant force of 0.002kBT/Rc is applied to each solvent particle in the X direction, and the simulation for several hundred thousand steps is performed until a steady flow is reached with minor fluctuations of velocities around the averaged velocity profile vs(z). Panels a−c in Figure 8 show the profiles of the solvent velocity along with PB and solvent densities across the channel cross-section at different solvent compositions. The computed solvent velocity vs(z) shown by dashed lines depends on the magnitude of the applied force. The normalized velocity, vn(z), reduced to the mean solvent velocity (shown by the solid lines) does not depend on the magnitude of the applied force (eq 7), as shown in the SI (section 2). The use of the normalized velocities is practical for comparison of the role of solvent composition on the fluid flow in the chromatographic columns at a given volumetric flow rate. The normalization of vn(z) implies that the mean reduced velocity ⟨vn(z)⟩ = 1, ensuring the same volumetric flux at different solvent compositions. It is worth noting that PB conformation at a given solvent composition and respective PB density profile are not affected by the flow (see SI section 2); this feature was also shown experimentally.65,66 At the same time, the fluid velocity is sensitive to the change of the solvent

⟨vsPA ⟩≈

2/3ρb vm(w − wPB) (ρb − ⟨ρPB ⟩)wPB + (w − wPB)ρb

(14)

from which we get back to eq 8 when ⟨ρPB⟩ = ρb (means no penetration of the solvent into the brush) and also when wPB = 0. When ⟨ρPB⟩ ≪ ρb, ⟨vsPA⟩ approaches the limiting value 2/3 vm(w − wPB)/w. Thus, as the mean velocity is smaller in the presence of PB, vn(z) becomes larger than 3/2. The mean solvent velocity of the PB-grated channel is determined approximately by the relative hydrodynamic width ((w − wPB)/w = wH/w) of the channel that determines the fraction of the solvent that flows. In Figure 8a, the normalized velocities become as high as close to 3 for the narrowest channel. Figure 8b presents the density and velocity distributions in a channel with PB of density 0.36 nm−2 for different solvent compositions. The hydrodynamic width of the channel in this case is increased for all xG because of the low grafting density. The magnitude of the solvent velocities also increases because of the lower PB density causing 1489

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir less resistance to the flow and the increased volume fraction of the solvent. On the other hand, vn(z) is lowered because of the larger hydrodynamic width of the channel. The use of a less repelling bad solvent B2 also causes changes in the velocity distributions (Figure 8c). In this case, systems 1 and 3 at xG = 0.99 are nearly the same because the bad solvent fraction is negligible, and as a result, the velocity and density profiles are similar. However, as xG decreases, the change of the PB conformation and the flow velocity are more moderate than in system 1. The hydrodynamic width of PB, wPB, is found by fitting vs(z) obtained from the simulations with parabolic approximation vsPA(z) at a given mean velocity

∫w

w

vsPA(z) dz =

PB

∫0

w

vs(z) dz

(15)

wPB can be iteratively found by solving eq 14 (see section 3 in the SI for details). The values of wPB obtained for various systems are shown in Table 2. Table 2. Hydrodynamic PB Thickness (wPB) at Different Good Solvent Fractions for Systems 1−3 wPB (Rc) xG

system 1

system 2

system 3

0.99 0.92 0.90 0.89 0.87 0.80 0.72 0.70 0.65 0.60

21.1 19.8 19.2

16.5 15.0

21.2

18.8 17.5

14.4 14.2 12.9

Figure 9. (a) Comparison of parabolic approximation (solid lines) with simulation results (circles, squares). Red, green, and blue represent solvent velocities vs(z) for xG = 80 in systems 1−3, respectively. Black, cyan, and magenta are the corresponding normalized profiles vn(z). (b) Normalized velocity profiles (solid lines) at xG = 0.99 in system 1 in channels of different widths w obtained using the parabolic approximation. Black, w = 33.4Rc; red, w = 50Rc; green, w = 100Rc; blue, w = 200Rc; magenta, w = 500Rc. The hydrodynamic thickness 21.12Rc is used. The solvent densities (dashes lines), which are taken to be constant in the bulk, are also shown.

20.1 19.5 19.4 19.0 18.3

of 3/2 that would be observed in the channels with solid walls without PB. 3.3. Henry Constant of NP Adhesion. We quantify NP− PB adhesion by the Henry constant calculated according to eq 6 from the free energy landscapes determined in the simulations with the GT method. The Henry constant is calculated by eq 6 as the excess of adsorption per unit substrate surface area using the reference bulk volume based on the PB hydrodynamic thickness wPB defined from the Poiseuille approximation of the solvent flow at a given solvent composition. The calculated Henry constant KH as a function of the good solvent fraction xG is presented in Figure 10. In the repulsive size exclusion regime, when the good solvent fraction is sufficiently high, the free energy is either positive or close to zero and the integrand in eq 10, e−A(z)/ kT ≤ 1. As a result, KH is negative. In a limiting case of very dense brush and purely repulsive interactions, KH ≈ −RNP, as in the HDC case, because the distances z < wPB + RNP from the substrate are effectively inaccessible to NP. In HDC,61 this size exclusion leads to shorter elution times for larger NPs: they are located closer to the channel center and move faster within the Poiseuille flow. Thus, we call the repulsive region of the NP− PB interaction at high good solvent fractions as the hydrodynamic, or size exclusion, mode, where KH is negative and decreases as the NP size increases. As the NP free energy profiles at lower good solvent fractions have regions of attraction with negative energy values, the

The hydrodynamic thickness wPB is the appropriate definition of the effective PB thickness when our purpose is to describe the flow through PB-grafted channels. It is worth noting that the hydrodynamic PB thickness is different from the effective geometrical thickness wG defined by approximating the PB density profile with the rectangular one assuming a uniform equilibrium PB density at z < wG, see section 3 of the SI for details. Figure 9a compares vs(z) and vn(z) obtained from the simulations with the corresponding parabolic approximations at xG = 0.80. Both solvent (circles) and normalized velocities (squares) showed very good agreement with their respective parabolic approximations for all systems (solid lines). Figure 9b presents the predicted velocity distributions using parabolic approximations for channels of different widths ranging from the width of the simulated channel w0 = 33.4Rc to 500Rc containing PB (system 1) at xG = 0.99. The solvent density ρs(z) in the channels wider than w0 was assumed equal to the computed solvent density in the channel of width w0 at z < w0 and set equal to the bulk density ρb at w0 < z < w. The solvent densities in Figure 9b are shown by dashed lines. The normalized velocities vnPA(z) are shown in Figure 9b by solid lines. Because they are constructed to provide the reduced flux equal to 1, the maximum normalized velocity in the channel center becomes smaller as the channel size increases. At w = 500Rc, the maximum value of vsPA(z) = 1.54, close to the value 1490

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

Figure 10. Solvent composition dependence of the Henry constant in systems (a) 1, (b) 2, and (c) 3. The Y axis represents the decimal logarithm of KH shifted by 10 nm to avoid negative values under logarithm in the size exclusion regime. Insets present KH−xG dependence to show the transition from the size exclusion regime with negative KH to the adsorption regime with positive KH.

Figure 11. Probability distributions P(z) of nanoparticles in the channel with PB grafting density 0.6 nm−2 (system 1) at good solvent fractions (a) 0.99, (b) 0.90, and (c) 0.80. The solid red, green, and blue curves represent P(z) of NPs with diameters of 16Rc, 12Rc, and 8Rc, respectively, and the dashed vertical lines correspondingly mark the mobility boundary wPB + RNP for each NP on PBs of both sides. The PB relative density ρPB′(z) (black) and the normalized velocity profiles vn(z) (magenta) are also shown.

surface chemistry can be achieved at least for select classes of NP mixtures. 3.4. Mean Velocity and Retention Time of NPs in PBGrafted Channels. To analyze NP separation in polymergrafted channels under flow, we follow the approximate approach adopted in HDC61,63,64 and assume that the NP velocity equals the solvent velocity at the position of the NP center of mass.17 Further, we employ the parabolic approximation of the flow velocity, introduce the hydrodynamic thickness of PB, wPB, and assume that NPs immersed (at least partially) in PB (z < wPB + RNP) are adsorbed and therefore immobile, whereas NPs located outside PB (that is, at z > wPB + RNP) move with the solvent velocity vs(z). As such, the mean NP velocity and retention time are calculated by convolution of the solvent velocity in parabolic approximation and the NP probability distribution according to eqs 12 and 13. In the absence of PB, this approximation reduces to the conventional HDC eqs 10 and 11. The probability density distribution P(z) of NP location in the channel is determined by the NP−PB free energy landscape for particular NP size and solvent composition according to eq 4. Figure 11a shows P(z) at xG = 0.99 in system 1 along with the PB relative density, ρPB′(z)= ρPB(z)/ρPBmax and normalized velocity vn(z) obtained from simulations. The assumed boundary of NP mobility, wPB + RNP, is shown as dashed lines. Because NPs are completely repelled from PB in such a good solvent, they tend to distribute uniformly in the bulk flow where the free energy is zero and the probability density is

exponential factor in the integrand of eq 10 becomes large and positive, and consequently, the Henry constants with the addition of bad solvent become positive and grow exponentially. This situation is referred as the adsorption regime of chromatography, in which the Henry constant KH is positive and increases with the NP size. As KH values vary with the good solvent fraction in an exponential manner, it is convenient to represent them in logarithmic scale. To avoid negative values under logarithm in the hydrodynamic regime, we shift KH by 10 nm and plot log(KH + 10 nm) against the good solvent fraction xG in Figure 10. The Henry constant values without such modification near the crossover are displayed in the insets to demonstrate the transition from the size exclusion regime with negative KH to the adsorption regime with positive KH. The transition from the size exclusion regime with negative KH in good solvents to the adsorption regime with positive KH at poorer solvent is observed in all systems modeled in this work. In each system, this transition happens in a relatively narrow range of solvent compositions: decreasing xG by approximately 0.02 reverses the order, the largest NPs that are most strongly repelled become most strongly adsorbed. From our results, we cannot pinpoint the exact value of xG where the Henry constant becomes NP size independent (that xG would correspond to the critical point of NP adsorption at PB). Nevertheless, the narrowness of the concentration range corresponding to the crossover in each particular system and relatively strong dependence of the location of the crossover region on the interactions indicate that NP separation by 1491

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

Figure 12. Transition from the size exclusion to adsorption regimes upon the decrease of solvent quality in systems (I) 1, (II) 2, and (III) 3. The channel half-width w = 70.43Rc = 50 nm. (a) Mean nanoparticle velocity ⟨vp⟩. (b) NP retention times reduced to solvent retention times τp/τs.. τp/τs < 1 in the size exclusion regime; τp/τs > 1 in the adsorption regime, and τp/τs ≈ 1 in the transition regime.

constant over z. It is clear from the plot that the NP probability density at z < wPB + RNP is nearly zero, whereas at z > wPB + RNP, it is constant. The hydrodynamic width of PB, wPB serves as a wall preventing penetration of NPs beyond z = wPB + RNP. In addition, the size exclusion constrains the larger particle toward the channel center, where it moves faster with the Poiseuille-type flow, whereas the smaller particles sample a wider space and wider range of velocities and, therefore, move slower on average. This picture is similar to that in HDC. Figure 11b shows probability distributions for NP location at xG = 0.90. Here, adsorption causes penetration of NPs in the interior of PB beyond the mobility boundary. However, the chances of finding the NP in the solvent bulk are still nonzero as the NP adsorption is weak. The overall velocity of NP in this case is determined by the fraction of NPs retained within PB where the NP velocity is zero. At xG = 0.80, as shown in Figure 11c, the NP probability density is concentrated within PB at z < wPB + RNP. Under these conditions, NPs are mainly retained in PB in the adsorbed state with a small probability to desorb. In this adsorption mode, the probability of desorption is determined by the depth of the adsorption well and decreases with the increase in NP size. Respectively, the mean velocity of NP decreases with their size. Figure 12a depicts the mean NP velocity ⟨vp⟩ calculated using eq 9 as a function of the solvent composition in systems (I) 1, (II) 2, and (III) 3. ⟨vp⟩ is calculated using the parabolic approximation assuming a channel haft-width, w = 70.43Rc, which is equivalent to a real channel width of 100 nm. In each case, ⟨vp⟩ decreases with the bad solvent fraction as the adsorption interactions gradually strengthen and attract NPs toward the PB interior. This effectively makes the probability of NP location in the mobile bulk fluid lower and therefore

restricts its mobility. In all three systems, the dependence of ⟨vp⟩ on the NP size undergoes a transition from size exclusion to adsorption modes. At higher xG, the mean velocities are larger for larger NPs due to the size exclusion of the NPs from the effectively repulsive PB as shown in Figure 11. At low xG, the NP probability density is determined by the adsorption energy, the larger particles with larger free energy of adsorption have a smaller probability to be located in the mobile solvent bulk at any given moment. Thus, the mean velocities are smaller for larger particles. When the adsorption is very strong, the NP velocity is nearly zero (Figure 12a I and II) but still decreases with the NP size. At the intermediate values of xG, the adsorption is weak, and the size dependence of ⟨vp⟩ arises from both size exclusion and adsorption effects, which approximately counter balance each other. Figure 12b shows the logarithmic values of NP retention times reduced by the solvent retention times (τp/τs), which is the reciprocal of ⟨vp⟩ (eq 13). In all cases, the retention times in the adsorption regime indicate that NPs are retained longer than the solvent (τp /τs > 1), whereas in the size exclusion regime, the solvent is retained longer (τp/τs < 1). In the transition regime, τp/τs ≈ 1 as shown in the insets in Figure 12b, the solvent and NPs elute within an approximately equal amount of time. The size exclusion mode of NP separation in the PB-grafted channel realized at high fractions of good solvent is similar to the conventional HDC with solid walls without PB. As mentioned in section 2.3, one may expect that NP mean velocity in the size exclusion mode approaches ⟨vpHDC⟩ given by eq 10 with δ = RNP/(w − wPB). In Figure 12a, the ⟨vp⟩ values for different NPs at xG = 0.99 in system 1 are 1.323, 1.387, and 1.450 for NPs of radius 4Rc, 6Rc, and 8Rc, respectively. Using the values wPB = 21.1Rc and w = 70.4Rc, and using δ = RNP/(w 1492

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

Figure 13. Ratio of retention times of NPs of different sizes to the solvent retention time as a function of the channel width in the (a, top) size exclusion and (b, bottom) adsorption regimes in systems (I) 1, (II) 2 and (III) 3. Colors: red, RNP = 8Rc; green, RNP = 6Rc; and blue, RNP = 4Rc.

− wPB) in eq 10, we obtain ⟨vpHDC⟩ of 1.079, 1.117, and 1.154 for the respective NPs. Thus, HDC values are ∼20% smaller. Similar difference exists between the mean velocity values in the size exclusion mode for systems 2 and 3 and the corresponding HDC values. The HDC values are smaller in all cases because, as described in Figure 9, the normalized velocities in the PBgrafted channels are larger than the velocities in the channels without PB due to the immobile solvent retained by the PB. ⟨vp⟩ approaches ⟨vpHDC⟩ as the channel size becomes larger. Figure 13 illustrates the effect of the channel width on the retention time in the size exclusion and adsorption regimes. In the size exclusion regime at xG = 0.99, the retention time as a function of the channel size is shown in Figure 13a. Here, NPs are not retained by the PB and are thus mobile, and the NP retention time increases as the channel width increases. This is consistent with the expectations in the HDC case, as according to eq 10, NP velocity decreases with the channel width due to the decrease of factor δ. In comparison with τpHDC (eq 11), τp is smaller because ⟨vp⟩ is larger than ⟨vpHDC⟩. Figure 13b presents the retention time in the adsorption mode. In this mode, τp increases with w in narrow channels in the case of NPs of radius 6Rc and 8Rc but decreases with channel width in wider channels. The range of w in which τp increases also depends on the NP size. For RNP = 8Rc, retention time increases with w until w reaches approximately 100Rc whereas for RNP = 6Rc, τp is found to increase only up to w ≈ 50Rc. In the case of the smallest NP of RNP= 4Rc, the similar behavior is therefore expected in smaller channels that are not shown in the figure. The increase in the retention time with the channel width, as in the size exclusion regime, is due to the decrease in the NP velocity as discussed above. However, the retention time decreases with w in larger channels. The reason for this behavior is the following: in the adsorption mode, NPs are

retained by the PB and stay immobile in the PB. The retention time in this case depends on the probability of NP desorption. The NP probability density in PB decreases as w increases (see eq 4), and therefore, the probability of desorption increases with channel width. Consequently, the presence of NPs in the bulk increases where they move with the flow and the retention time decreases.

4. CONCLUSIONS This work explores the effects of solvent quality on NP adhesion to polymer brushes and NP flow through polymergrafted channels. These physical phenomena are important for multiple applications and, in particular, to the problem of chromatographic separation and purification of chemically functionalized nanoparticles that is the practical focus of this work. We show that the energy of NP−PB adhesion and, respectively, NP partition between the stationary phase represented by PB-grafted pore surface and the mobile phase of flowing solvent critically depend on the solvent quality with respect to the polymer. To explicitly control the solvent quality, we employ a binary mixture of good and bad solvent components and vary its composition. Using coarse-grained DPD simulations, we demonstrate the characteristic features of NP−PB interactions in binary solvent drawing on an example of a model system parametrized to mimic polyisoprene natural rubber brush in the mixture of benzene (good solvent) and acetone (bad solvent) and NPs functionalized by short chains. NP−PB adhesion is governed by adsorption of polymer chains to the NP surface and entropic repulsion caused by polymer chain confinement between NP and the channel wall. Both factors depend on the solvent composition, variation of which causes (a) contraction or expansion of PB and (b) alteration of the NP affinity to the polymer. At good solvent 1493

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

NP moves with the solvent velocity. The mean NP velocity is calculated by the convolution of the solvent velocity profile with the probability density of NP location within the channel; the latter is determined by the Boltzmann distribution within the free energy landscape independently calculated using the GT method. Thus, calculated mean fluid velocities allow us to compare the respective retention times depending on the solvent composition provided the volumetric flow rate is fixed to mimic chromatographic experiments. The obtained dependence convincingly confirms the transition from the size exclusion to adsorption regimes and respective reversal of the sequence of elution occurring with the decrease of solvent quality in a narrow range of solvent compositions. The observed reversal of the sequence of elution is reminiscent of the critical condition in polymer chromatography at which the retention time is molecular weight independent. In polymer chromatography, this transition occurs at the so-called critical point of adsorption (CPA), and the corresponding separation regime, often called liquid chromatography at critical conditions (LCCC), has been widely applied to size-independent separation polymers by their chemistry and architecture in the isocratic and gradient elution modes.25−27 Although we cannot not identify an analogue of CPA which is characteristic for polymer adsorption, the transition from size-exclusion to adsorption regimes of NP separation is found in a narrow interval of solvent compositions. This finding suggests the possibility of identifying the separation conditions (composition of solvent that depends on the NP chemistry) at which NPs with like surface properties elute together regardless of their size. The latter has important practical implications: NPs can be separated by surface chemistry rather than by their size, employing the gradient mode of elution with controlled variation of solvent composition, similarly to that achieved in LCCC of polymers. This regime, which might be named nanoparticle chromatography at critical conditions (NPCCC), is highly desirable for characterization, purification, and quality control of chemically modified nanoparticles. Further modeling and experimental work should be performed to identify the dependence of the critical conditions of NP separation on NP functionalization in given polymer−solvent systems. It is worth noting that we considered a limited range of nanoparticle sizes (the diameter of the solid particle core ranged from 5.5 to 11.1 nm) that was sufficient to demonstrate the transition between the size exclusion and adsorption regimes of nanoparticle separation that is the main result of this paper. It would be very desirable to extend the range of particle sizes but such simulations were prohibitively expensive given the scope of this work. At a limit of large particles (that is, when particle diameter in substantially larger than the brush thickness), the adhesion can be quantified in terms of the disjoining pressure between the particle surface and the substrate coated by the polymer brush. This makes an interesting point for further investigations.

conditions, PB is expanded and NP is repelled entropically from PB with the force that increase with the NP size. This sizeexclusion regime is similar to that of hydrodynamic chromatography (HDC) on capillary or porous columns: the retention time reduces with the particle size so that the larger particles elute faster. This trend dramatically changes as the solvent quality worsens upon decrease of the good solvent fraction causing contraction and densification of PB. Therefore, the polymer adsorption to NP becomes stronger due to the densification of PB and expulsion of the bad solvent component from the NP−PB contact region. This facilitates the enthalpic attraction of NP to PB, which at sufficiently poor solvent conditions prevails over the entropic repulsion. This adsorption regime is similar to that in liquid chromatography of polymers: enthalpic attraction increases with the particle size, causing larger retention time for larger particles. In our simulations, we demonstrate that, with the increase in the bad solvent fraction, NP separation exhibits a transition from the size exclusion regime with larger NPs having shorter retention time to the adsorption regime with smaller NPs having shorter retention time. Using extensive dissipative particle dynamics simulations in conjunction with the ghost tweezers free energy calculation technique, we quantified the specifics of NP−PB adhesion interactions by calculating the free energy landscapes of functionalized NPs within PB-grafted channels depending on the solvent composition for three characteristic systems differing by PB grafting density and solvent−polymer affinity. The free energy landscape determines the probability of NP location at given distance to the surface, position of equilibrium adhesion state with partial immersion of NP into PB, and the Henry constant that characterizes adsorption equilibrium. We find that, for all systems considered, the Henry constant is negative and decreases with the NP size at good solvent conditions that reflects the predominance of entropic repulsion characteristics to the hydrodynamic size exclusion regime. Upon the increase of the bad solvent fraction, the Henry constant increases and becomes positive, indicating the transition into the adsorption regime with larger particles more strongly adhered to PB. This transition for all particle sizes occurs in a narrow interval of solvent compositions so that the transition composition can be estimated within a reasonable error. To analyze the effects of solvent quality on fluid flow through polymer-grafted pores, we simulated the solvent flow in slitshaped channels at various compositions of good and bad solvent. Addition of the bad solvent leads to PB contraction and a respective increase of the channel volume available for the flow. We find that, for all systems considered, the solvent flow can be approximated by a Poiseuillian parabolic profile with zero velocity inside PB within a certain distance from the solid substrate that was called the hydrodynamic width of PB, wPB. As such, PB acts as a “solvated” wall of width wPB beyond which the velocity profile is similar to that in a channel with solid walls. We considered channels of different widths and concluded that the PB hydrodynamic width does not depend on the width of the channel and mean fluid velocity. This conclusion allows us to predict the solvent flow in wider channels from the simulation performed in smaller channels, which require less computational time. The solvent velocity profile in the parabolic approximation is further used to estimate the mean NP velocity. To this end, we assume, as is conventionally done in hydrodynamic chromatography, that

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.langmuir.7b02914. Brief sketches of the effects of the presence of nanoparticles and the magnitude of the applied forces 1494

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir

■

(14) Prokhorova, S. A.; Kopyshev, A.; Ramakrishnan, A.; Zhang, H.; Rühe, J. Can polymer brushes induce motion of nano-objects? Nanotechnology 2003, 14 (10), 1098. (15) Arita, T.; Yoshimura, T.; Adschiri, T. Size exclusion chromatography of quantum dots by utilizing nanoparticle repelling surface of concentrated polymer brush. Nanoscale 2010, 2 (8), 1467− 1473. (16) Wei, G.-T.; Liu, F.-K. Separation of nanometer gold particles by size exclusion chromatography. J. Chromatogr. A 1999, 836 (2), 253− 260. (17) Striegel, A. M.; Brewer, A. K. Hydrodynamic Chromatography. Annu. Rev. Anal. Chem. 2012, 5 (1), 15−34. (18) Malysheva, A.; Lombi, E.; Voelcker, N. H. Bridging the divide between human and environmental nanotoxicology. Nat. Nanotechnol. 2015, 10 (10), 835−844. (19) Kowalczyk, B.; Lagzi, I.; Grzybowski, B. A. Nanoseparations: Strategies for size and/or shape-selective purification of nanoparticles. Curr. Opin. Colloid Interface Sci. 2011, 16 (2), 135−148. (20) Lee, J.-S.; Stoeva, S. I.; Mirkin, C. A. DNA-Induced SizeSelective Separation of Mixtures of Gold Nanoparticles. J. Am. Chem. Soc. 2006, 128 (27), 8899−8903. (21) Albanese, A.; Tang, P. S.; Chan, W. C. W. The Effect of Nanoparticle Size, Shape, and Surface Chemistry on Biological Systems. Annu. Rev. Biomed. Eng. 2012, 14 (1), 1−16. (22) Kim, S. T.; Saha, K.; Kim, C.; Rotello, V. M. The Role of Surface Functionality in Determining Nanoparticle Cytotoxicity. Acc. Chem. Res. 2013, 46 (3), 681−691. (23) Paget, V.; Dekali, S.; Kortulewski, T.; Grall, R.; Gamez, C.; Blazy, K.; Aguerre-Chariol, O.; Chevillard, S.; Braun, A.; Rat, P.; Lacroix, G. Specific Uptake and Genotoxicity Induced by Polystyrene Nanobeads with Distinct Surface Chemistry on Human Lung Epithelial Cells and Macrophages. PLoS One 2015, 10 (4), e0123297. (24) Guttman, C. M.; Di Marzio, E. A.; Douglas, J. F. Influence of Polymer Architecture and Polymer−Surface Interaction on the Elution Chromatography of Macromolecules through a Microporous Media. Macromolecules 1996, 29 (17), 5723−5733. (25) Brun, Y. The Mechanism of Copolymer Retention in Interactive Polymer Chromatography. I. Critical Point of Adsorption for Statistical Copolymers. J. Liq. Chromatogr. Relat. Technol. 1999, 22 (20), 3027−3065. (26) Brun, Y. The Mechanism of Copolymer Retention in Interactive Polymer Chromatography. II. Gradient Separation. J. Liq. Chromatogr. Relat. Technol. 1999, 22 (20), 3067−3090. (27) Brun, Y.; Alden, P. Gradient separation of polymers at critical point of adsorption. J. Chromatogr. A 2002, 966 (1−2), 25−40. (28) Marko, J. F. Polymer brush in contact with a mixture of solvents. Macromolecules 1993, 26 (2), 313−319. (29) Cimino, R. T.; Rasmussen, C. J.; Brun, Y.; Neimark, A. V. Critical conditions of polymer adsorption and chromatography on non-porous substrates. J. Colloid Interface Sci. 2016, 474, 25−33. (30) Cimino, R. T.; Rasmussen, C. J.; Brun, Y.; Neimark, A. V. Mechanisms of chain adsorption on porous substrates and critical conditions of polymer chromatography. J. Colloid Interface Sci. 2016, 481, 181−193. (31) Halperin, A.; Kröger, M.; Zhulina, E. B. Colloid-Brush Interactions: The Effect of Solvent Quality. Macromolecules 2011, 44 (9), 3622−3638. (32) Zhang, X.; Chen, Y.; Qu, L.; Yan, D. Effects of attractive colloids on the phase separation behaviors of binary polymer blends. J. Chem. Phys. 2013, 139 (7), 074902. (33) Lewis, T.; Ganesan, V. Interactions between Grafted Cationic Dendrimers and Anionic Bilayer Membranes. J. Phys. Chem. B 2013, 117 (33), 9806−9820. (34) LoVerso, F.; Egorov, S. A.; Binder, K. Interaction Between Polymer Brush-Coated Spherical Nanoparticles: Effect of Solvent Quality. Macromolecules 2012, 45 (21), 8892−8902. (35) Surve, M.; Pryamitsyn, V.; Ganesan, V. Polymer-bridged gels of nanoparticles in solutions of adsorbing polymers. J. Chem. Phys. 2006, 125 (6), 064903.

on the solvent velocity profiles, iteration procedure for estimating the hydrodynamic PB thickness, and definition and numerical values of the geometrical thickness of PBs (PDF)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Alexander V. Neimark: 0000-0001-6817-737X Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was funded by the National Science Foundation CBET grant “GOALI: Theoretical Foundations of Interaction Nanoparticle Chromatography”, No. 510993, and used computational resources of Extreme Science and Engineering Discovery Environment (XSEDE) supported by the National Science Foundation Grant ACI-1548562.

■

REFERENCES

(1) Peer, D.; Karp, J. M.; Hong, S.; Farokhzad, O. C.; Margalit, R.; Langer, R. Nanocarriers as an emerging platform for cancer therapy. Nat. Nanotechnol. 2007, 2 (12), 751−760. (2) Roy, I.; Ohulchanskyy, T. Y.; Pudavar, H. E.; Bergey, E. J.; Oseroff, A. R.; Morgan, J.; Dougherty, T. J.; Prasad, P. N. CeramicBased Nanoparticles Entrapping Water-Insoluble Photosensitizing Anticancer Drugs: A Novel Drug−Carrier System for Photodynamic Therapy. J. Am. Chem. Soc. 2003, 125 (26), 7860−7865. (3) Tanaka, K.; Chujo, Y. Design of functionalized nanoparticles for the applications in nanobiotechnology. Adv. Powder Technol. 2014, 25 (1), 101−113. (4) Anker, J. N.; Hall, W. P.; Lyandres, O.; Shah, N. C.; Zhao, J.; Van Duyne, R. P. Biosensing with plasmonic nanosensors. Nat. Mater. 2008, 7 (6), 442−453. (5) Klajn, R.; Browne, K. P.; Soh, S.; Grzybowski, B. A. Nanoparticles That “Remember” Temperature. Small 2010, 6 (13), 1385−1387. (6) Nakanishi, H.; Grzybowski, B. A. Supercapacitors Based on Metal Electrodes Prepared from Nanoparticle Mixtures at Room Temperature. J. Phys. Chem. Lett. 2010, 1 (9), 1428−1431. (7) Bliznyuk, V.; Ruhstaller, B.; Brock, P. J.; Scherf, U.; Carter, S. A. Self-Assembled Nanocomposite Polymer Light-Emitting Diodes with Improved Efficiency and Luminance. Adv. Mater. 1999, 11 (15), 1257−1261. (8) Park, J. H.; Lim, Y. T.; Park, O. O.; Kim, J. K.; Yu, J.-W.; Kim, Y. C. Polymer/Gold Nanoparticle Nanocomposite Light-Emitting Diodes: Enhancement of Electroluminescence Stability and Quantum Efficiency of Blue-Light-Emitting Polymers. Chem. Mater. 2004, 16 (4), 688−692. (9) Hasobe, T.; Imahori, H.; Kamat, P. V.; Ahn, T. K.; Kim, S. K.; Kim, D.; Fujimoto, A.; Hirakawa, T.; Fukuzumi, S. Photovoltaic Cells Using Composite Nanoclusters of Porphyrins and Fullerenes with Gold Nanoparticles. J. Am. Chem. Soc. 2005, 127 (4), 1216−1228. (10) Sun, B.; Marx, E.; Greenham, N. C. Photovoltaic Devices Using Blends of Branched CdSe Nanoparticles and Conjugated Polymers. Nano Lett. 2003, 3 (7), 961−963. (11) Tseng, R. J.; Tsai, C.; Ma, L.; Ouyang, J.; Ozkan, C. S.; Yang, Y. Digital memory device based on tobacco mosaic virus conjugated with nanoparticles. Nat. Nanotechnol. 2006, 1 (1), 72−77. (12) Ayres, N. Polymer brushes: Applications in biomaterials and nanotechnology. Polym. Chem. 2010, 1 (6), 769−777. (13) Kreer, T. Polymer-brush lubrication: a review of recent theoretical advances. Soft Matter 2016, 12 (15), 3479−3501. 1495

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496

Article

Langmuir (36) Chen, K.; Ma, Y.-q. Interactions between Colloidal Particles Induced by Polymer Brushes Grafted onto the Substrate. J. Phys. Chem. B 2005, 109 (37), 17617−17622. (37) Egorov, S. A.; Binder, K. Effect of solvent quality on the dispersibility of polymer-grafted spherical nanoparticles in polymer solutions. J. Chem. Phys. 2012, 137 (9), 094901. (38) Guskova, O. A.; Pal, S.; Seidel, C. Organization of nanoparticles at the polymer brush-solvent interface. Europhys. Lett. 2009, 88 (3), 38006. (39) Kim, J. U.; Matsen, M. W. Repulsion Exerted on a Spherical Particle by a Polymer Brush. Macromolecules 2008, 41 (1), 246−252. (40) Zhang, Q.; Xiang, X. Adsorption of a spherical nanoparticle in polymer brushes: Brownian dynamics investigation. Phys. A 2013, 392 (18), 3857−3862. (41) Curk, T.; Martinez-Veracoechea, F. J.; Frenkel, D.; Dobnikar, J. Collective ordering of colloids in grafted polymer layers. Soft Matter 2013, 9 (23), 5565−5571. (42) Romeis, D.; Sommer, J.-U. Conformational switching of modified guest chains in polymer brushes. J. Chem. Phys. 2013, 139 (4), 044910. (43) Vogiatzis, G. G.; Theodorou, D. N. Structure of Polymer Layers Grafted to Nanoparticles in Silica−Polystyrene Nanocomposites. Macromolecules 2013, 46 (11), 4670−4683. (44) Murat, M.; Grest, G. S. Molecular Dynamics Simulations of the Force between a Polymer Brush and an AFM Tip. Macromolecules 1996, 29 (25), 8282−8284. (45) Milchev, A.; Dimitrov, D. I.; Binder, K. Excess free energy of nanoparticles in a polymer brush. Polymer 2008, 49 (17), 3611−3618. (46) Bonet Avalos, J.; Johner, A.; Dıez-Orrite, S. Analysis of polymer adsorption onto colloidal particles. Eur. Phys. J. E: Soft Matter Biol. Phys. 2006, 21 (4), 305−317. (47) Yang, S.; Yan, D. D.; Shi, A. C. Structure of adsorbed polymers on a colloid particle. Macromolecules 2006, 39 (12), 4168−4174. (48) Merlitz, H.; Wu, C. X.; Sommer, J. U. Inclusion Free Energy of Nanoparticles in Polymer Brushes. Macromolecules 2012, 45 (20), 8494−8501. (49) Cheng, J.; Vishnyakov, A.; Neimark, A. V. Morphological Transformations in Polymer Brushes in Binary Mixtures: DPD Study. Langmuir 2014, 30 (43), 12932−12940. (50) Hua, Y.; Zhang, D.; Zhang, L. Compression-driven migration of nanoparticles in semiflexible polymer brushes. Polymer 2016, 83, 67− 76. (51) Nasrabad, A. E.; Jasnow, D.; Zilman, A.; Coalson, R. D. Precise control of polymer coated nanopores by nanoparticle additives: Insights from computational modeling. J. Chem. Phys. 2016, 145 (6), 064901. (52) Cao, Q. Q.; Zuo, C. C.; Li, L. J.; Li, Y. J.; Yang, Y. Translocation of nanoparticles through a polymer brush-modified nanochannel. Biomicrofluidics 2012, 6 (3), 10. (53) Cheng, J.; Vishnyakov, A.; Neimark, A. V. Adhesion of nanoparticles to polymer brushes studied with the ghost tweezers method. J. Chem. Phys. 2015, 142 (3), 034705. (54) Español, P.; Warren, P. Statistical Mechanics of Dissipative Particle Dynamics. Europhys. Lett. 1995, 30 (4), 191. (55) 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. (56) Vishnyakov, A.; Talaga, D. S.; Neimark, A. V. DPD Simulation of Protein Conformations: From α-Helices to β-Structures. J. Phys. Chem. Lett. 2012, 3 (21), 3081−3087. (57) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117 (1), 1−19. (58) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual molecular dynamics. J. Mol. Graphics 1996, 14 (1), 33−38. (59) Matthews, J. N. A. Commercial optical traps emerge from biophysics labs. Phys. Today 2009, 62 (2), 26−28. (60) Yang, S.; Neimark, A. V. Critical conditions of polymer chromatography: An insight from SCFT modeling. J. Chem. Phys. 2013, 138 (24), 244903.

(61) DiMarzio, E. A.; Guttman, C. M. Separation by Flow. Macromolecules 1970, 3 (2), 131−146. (62) Small, H. Hydrodynamic chromatography a technique for size analysis of colloidal particles. J. Colloid Interface Sci. 1974, 48 (1), 147− 161. (63) Brenner, H.; Gaydos, L. J. The constrained brownian movement of spherical particles in cylindrical pores of comparable radius. J. Colloid Interface Sci. 1977, 58 (2), 312−356. (64) Prieve, D. C.; Hoysan, P. M. Role of colloidal forces in hydrodynamic chromatography. J. Colloid Interface Sci. 1978, 64 (2), 201−213. (65) Baker, S. M.; Smith, G. S.; Anastassopoulos, D. L.; Toprakcioglu, C.; Vradis, A. A.; Bucknall, D. G. Structure of Polymer Brushes under Shear Flow in a Good Solvent. Macromolecules 2000, 33 (4), 1120−1122. (66) Ivkov, R.; Butler, P. D.; Satija, S. K.; Fetters, L. J. Effect of Solvent Flow on a Polymer Brush: A Neutron Reflectivity Study of the Brush Height and Chain Density Profile. Langmuir 2001, 17 (10), 2999−3005.

1496

DOI: 10.1021/acs.langmuir.7b02914 Langmuir 2018, 34, 1481−1496