In polymer nanoparticle composites (PNCs) with attractive interactions between nanoparticles (NPs) and polymers, a bound layer of the polymer forms on the NP surface, with significant effects on the macroscopic properties of the PNCs. The adsorption and wetting behaviors of polymer solutions in the presence of a solid surface are critical to the fabrication process of PNCs. In this study, we use both classical density functional theory (cDFT) and molecular dynamics (MD) simulations to study dilute and semi-dilute solutions of short polymer chains near a solid surface. Using cDFT, we calculate the equilibrium properties of polymer solutions near a flat surface while varying the solvent quality, surface–fluid interactions, and the polymer chain lengths to investigate their effects on the polymer adsorption and wetting transitions. Using MD simulations, we simulate polymer solutions near solid surfaces with three different curvatures (a flat surface and NPs with two radii) to study the static conformation of the polymer bound layer near the surface and the dynamic chain adsorption process. We find that the bulk polymer concentration at which the wetting transition in the poor solvent system occurs is not affected by the difference in surface–fluid interactions; however, a threshold value of surface–fluid interaction is needed to observe the wetting transition. We also find that with good solvent, increasing the chain length or the difference in the surface–polymer interaction relative to the surface–solvent interaction increases the surface coverage of polymer segments and independent chains for all surface curvatures. Finally, we demonstrate that the polymer segmental adsorption times are heavily influenced only by the surface–fluid interactions, although polymers desorb more quickly from highly curved surfaces.

## I. INTRODUCTION

Polymer nanoparticle composites (PNCs) exhibit both macroscopic and microscopic properties that are distinct from those of pure polymers.^{1–17} In PNCs with attractive interactions between the nanoparticle (NP) surfaces and the polymers, polymers adsorb to the surface and form a bound layer. This polymer bound layer has been shown to be responsible for changes in the properties and performance of the PNCs and has been used to prevent NP aggregation^{18,19} and to improve mechanical^{20,21} and transport properties of PNCs.^{22–24} Significant progress has been made in understanding the multi-scale structure and segmental dynamics of the bound layer polymers in PNCs.^{4,9,23,25–28} First, the bound layer thickness in solution and in melt has been measured using various experimental techniques.^{4,25,26} The amount of excess adsorption of polymers is affected by the strength of the interaction between the surface and the polymer^{29} and by the bulk polymer concentration.^{30} The density of polymers in the bound layer is impacted by the NP–polymer interaction strength^{12,31,32} and the molecular weight of the polymers.^{12,33,34} Finally, the surface curvature plays a significant role in the conformation of polymers in the bound layer.^{5,7,13,16,27,35–40} In these experimental studies, the polymer chains have been shown to be compressed perpendicularly to the surface,^{36} which is in qualitative agreement with theoretical predictions^{41,42} and simulation results.^{16,27,36,43}

A common fabrication process of PNCs involves dispersing NPs in a polymer solution prior to annealing, during which the behavior of polymer chains near the NPs can influence the resulting structure and properties of the dried PNCs.^{18} Therefore, the solvents used in the fabrication process can result in significant differences in macroscopic properties from similar samples produced in different laboratories.^{18,44} It is widely accepted that in equilibrium, the bound layer of polymers on a surface in a polymer solution is not immobile, but instead, there is a continuous process of polymer chain exchanges between the bound layer and the bulk solution. Nevertheless, there is still some disagreement in the literature as to whether the bound layer formation is reversible.^{45–48} Thus, understanding the chain adsorption and exchange process is critical to understanding the various properties of the bound layer. Previous research on the kinetics of polymer chain adsorption has identified three separate stages in this process.^{49–51} First, the polymer chain diffuses close to the adsorbing surface, and when the chain approaches the surface, its conformation becomes distorted.^{51} Next, the polymer chain attaches itself onto the surface, then the segments rearrange, and the polymer flattens, which lowers the potential energy of the system, and can displace previously adsorbed polymers. This step of the process is slow and history-dependent, during which loops and trains can form on the surface.^{49,50,52} Finally, chain detachment occurs at a time scale exponentially proportional to the chain length.^{49} When the interaction between the surface and the polymer is favorable, the formation of a bound layer of polymers on the surface is observed, and the strength of the interaction affects the amount of the adsorbed polymer^{18,52,53} and the chain diffusivity near the surface.^{39}

A fundamental understanding of the static and dynamic properties of bound layer polymers remains elusive, especially at the molecular level. Most of the previous work focuses on adsorption of polymer chains onto smooth surfaces; in experimental systems of PNCs, the roughness of the solid surface cannot be ignored, especially in the case of small NPs with large surface area to volume ratios.^{54,55} In addition, the temporal properties of the bound layer, such as the time scales of its formation and the chain exchange processes, have not been well established.^{14,32,35,46} In this paper, we study the wetting and adsorption–desorption transitions in polymer solutions in contact with a solid surface using classical density functional theory (cDFT) for fluids and investigate polymer chain-scale adsorption dynamics in the bound layer using molecular dynamics (MD) simulations. A schematic of the systems we investigate is shown in Fig. 1. A common quantity used to describe the size of the polymer chain relative to the size of the nanoparticle is *R*_{g}/*r*_{NP}, where *R*_{g} is the ideal radius of gyration of the polymer and *r*_{NP} is the nanoparticle radius. We vary *R*_{g}/*r*_{NP} by varying the curvature of the solid surfaces and the polymer chain lengths, and we also vary solvent quality and the strength of polymer–surface interaction to examine their effects on the static and dynamic properties of the polymer bound layer. We use cDFT to construct phase diagrams of bulk polymer solutions with two polymer chain lengths, *N* = 10 and *N* = 40. Using two solvent qualities that represent good and poor solvents, we investigate the adsorption and wetting behavior, respectively. Combining with coarse-grained MD simulations, we explore the polymer density profiles and the adsorption dynamics as a function of the curvature, chain length, and polymer–surface interaction strength. The results from our study provide guidelines for experimental designs to probe bound layer conformations and dynamics in polymer solutions with NPs.

## II. METHODS

For both cDFT calculations and MD simulations, we use dilute to semi-dilute polymer solutions with low loading of NPs to investigate polymer adsorption in the absence of confinement and interactions with multiple surfaces. By varying the solvent quality (*ɛ*_{ps}) and the surface–fluid interactions, we mimic experimental systems containing NPs with different surface properties. For the fluid–fluid interactions, we use the 12-6 Lennard-Jones (LJ) potential $uf1,f2(r)=4\epsilon f1,f2\sigma f1,f2r12\u2212\sigma f1,f2r6$, where *f*_{1} and *f*_{2} indicate the two types of fluid particles in our system. We use the same polymer–polymer, solvent–solvent self-interaction strength, *ɛ*_{pp} = *ɛ*_{ss} = 1.0 kT, in both cDFT calculations and MD simulations. For the attractive surface–fluid interactions in our cDFT calculations, we use the 9-3 LJ potential $uwf(r)=332\epsilon wf\sigma wfr9\u2212\sigma wfr3$, where *w* and *f* indicate the surface and fluid (either polymer or solvent) particles. The surface–fluid interactions in our MD simulations are modeled using the 12-6 LJ potential. To ensure that the NPs are miscible in solution, we choose the surface–solvent interactions to be attractive throughout our study. We fix the surface–solvent interaction strength to be *ɛ*_{ws} = 1.25 kT for the cDFT calculations. For the MD simulations, we fix *ɛ*_{ws} = 0.65 kT to match the depth of the potential well between MD simulations and cDFT calculations. The details of this procedure are described in Sec. II B. All non-bonded interactions in both the cDFT calculations and the MD simulations have a cutoff length of *r*_{cut} = 2.5*σ*. Additionally, we define Δ*ɛ*_{wf} to be the difference in surface–fluid interactions as Δ*ɛ*_{wf} = *ɛ*_{wp} − *ɛ*_{ws}, where *ɛ*_{wp} is the surface–polymer interaction, so that Δ*ɛ*_{wf} > 0 represents a more favorable surface–polymer interaction compared to the surface–solvent interaction. In both cDFT calculations and MD simulations, lengths are in units of *σ* and energies are in units of *kT*. Note that the densities are in units of the number of segments per volume, *σ*^{−3}. All polymer monomers and solvent particles have a diameter of *σ* = 1.0. The units are in reduced LJ units: temperature *T* = *kT**/*ɛ* and time $\tau LJ=t*\epsilon /m\sigma 2$, where *m* represents the mass of a single LJ interaction site and *T** and *t** represent temperature and time measured in laboratory units.

We use two polymer chain lengths, *N* = 10 and *N* = 40, corresponding to ideal radii of gyration of *R*_{g} = 1.3*σ* and 2.6*σ*, respectively. Both chain lengths are below the entanglement chain length, *N*_{e}. For the cDFT calculations, we use only the flat surface; for the MD simulations, we use the flat surface, a NP with *r*_{NP} = 25*σ* and a NP with *r*_{NP} = 6*σ*. These combinations of polymer *R*_{g} and surface curvatures correspond to a range of *R*_{g}/*r*_{NP}, tabulated in Table I. Note that the overlap volume concentrations for these chain lengths in good solvent are approximately *ϕ** ∼ 0.16 and *ϕ** ∼ 0.05 for *N* = 10 and *N* = 40, respectively.

### A. Theoretical calculations

We use cDFT for fluids with bonded interactions and pairwise attractions to calculate the equilibrium thermodynamic properties of our system consisting of a solid surface in a polymer solution. The fluid-cDFT approach described here has previously been used to study polymers near surfaces.^{56–59} The theoretical approach is to find an equilibrium density profile for each fluid component $\rho \alpha r$ that minimizes the grand free energy functional $\Omega \rho \alpha r$,

where $F\rho \alpha r$ is the intrinsic Helmholtz free energy, *μ*_{α} is the chemical potential of species *α*, and $V\alpha r$ is an external field. The intrinsic Helmholtz free energy is expressed as

where the right-hand side terms are the intrinsic Helmholtz free energies for the ideal gas {$Fid\rho \alpha r$}, hard sphere {$Fhs\rho \alpha r$}, attractive interactions {$Fatt\rho \alpha r$}, and chain constraints {$Fch\rho \alpha r$}. The $Fhs\rho \alpha r$ term describes a reference fluid consisting of a simple fluid of hard spheres. Specifically, this model employs the White Bear fundamental measure theory (FMT) for the hard sphere energy term.^{60} The attractive interactions between particles and the bonding constraints are included as perturbations to the reference fluid. Finally, this fluid-cDFT model expresses the $Fch\rho \alpha r$ term using the modified interfacial statistical associating fluid theory (modified-i-SAFT) to describe chain constraints due to the presence of polymer bonds.^{61} The modified-i-SAFT model is based on an extension of Wertheim’s first order thermodynamic perturbation theory and is generally applicable to mixtures of inhomogeneous complex fluids.^{62}

Solving for the equilibrium density profiles for each fluid component that minimizes the grand free energy functional results in a system of nonlinear integral equations. In our study, we use the open-source software developed at Sandia National Laboratories, Tramonto, for the calculations.^{63} We also rely on a pseudo-arc-length continuation (ALC) algorithm developed previously to vary specific parameters and find coexisting solutions.^{64,65} A detailed description of the numerical treatment can be found in Appendixes B and C in the work of Frink *et al.*^{58,66} When varying the composition in a binary mixture, we need to keep the pressure constant to remain physically realistic. In modified-i-SAFT, the bulk pressure is an analytic function of polymer and solvent bulk densities, $p=f(\rho pb,\rho sb)$. This function can be inverted numerically to obtain the solvent bulk density as a function of polymer bulk density $\rho sb=g(\rho pb;p)$ at a constant pressure *p*.^{57} We find that we can fit the function $g(\rho pb;p)$ with a fourth-order polynomial (see the supplementary material), which is used to vary solvent density as a function of polymer density to keep the pressure constant.

We perform the cDFT calculations in 1D only due to the large computational cost required for 3D calculations, particularly, since our system includes attractive interactions. Our calculation box size is 80*σ*, with a mesh size of 0.05*σ*. One end of the computational domain has a semi-infinite surface to simulate a flat surface, and the other end has a reflective boundary condition, which means that the densities beyond the last computational node are set to the same densities an equal distance away from the boundary, i.e., for the *k*th node beyond the boundary, *ρ*_{n+k} = *ρ*_{n−k}, where *n* is the total number of nodes.

### B. Molecular dynamics simulation

In addition to the 1D cDFT calculations, we use MD simulations to extend our study to 3D to investigate the effects of surface curvature. We use the simulation package LAMMPS with periodic boundary conditions in three dimensions to simulate infinitely extended systems.^{67} All bonds are simulated with a harmonic bonding potential with a spring constant $Uijb=(kh/2)(rij\u2212\sigma )2$, where *k*_{h} = 2000 *ɛ*/*σ*^{2}, and the non-bonded interactions are the same as previously described. The above model for simulating fully flexible polymer chains and surface–fluid interaction has been well established.^{68} We use explicit particles with LJ interactions for the solvent. The parameters in the LJ potentials are set as the same as in the cDFT calculations: *σ*_{p} = *σ*_{s} = 1.0 and *ɛ*_{pp} = *ɛ*_{ps} = *ɛ*_{ss} = 1.0. We use a simulation time step size of *δt* = 0.002 *τ*_{LJ}.

To simulate the solid surfaces, we cut a rigid slab or sphere with defined width or radius from an amorphous LJ surface with a bulk density of *ρ* = 1.0. This process allows for the inclusion of roughness on the surface. The rigid slab simulating a flat surface has a thickness of 8*σ*, while the NPs have two radii, *r*_{NP} = 6*σ* and *r*_{NP} = 25*σ*. The centers of mass of the solid surfaces are fixed at the center of the simulation box, but they are allowed to rotate. Each particle in the surfaces can interact with the solvent and the polymer particles, which makes it necessary to tune the strength of the surface–fluid interaction to match the depth of the potential well with that of the 9-3 LJ potential used in the cDFT calculations. The resulting surface–solvent interaction strength is *ɛ*_{ws} = 0.65, and the surface–polymer interaction strength is varied to vary Δ*ɛ*_{wf}. The polymer and solvent densities are $\rho pb=0.1$ and $\rho sb=0.6$, respectively, and are calculated as $\rho \alpha b=n\alpha /(Vbox\u2212Vsurface)$, where *n*_{α} denotes the number of particles of a specific species and *V*_{box} and *V*_{surface} are the total simulation box volume and the surface volume, respectively. The polymer density corresponds to a polymer volume fraction of *ϕ* = 0.167, which is above the overlap concentration for both chain lengths, so the MD simulations are for semi-dilute solutions. Snapshots of the three systems are shown in Fig. 2.

Three independent configurations for each type of surface are equilibrated to provide estimates of the uncertainty in our calculations. Equilibration is confirmed with the observation of polymers reaching the diffusive regime via the calculation of the mean-squared displacement. After equilibration, particle positions are collected at multiple times during the simulation and used in our analyses.

## III. RESULTS

### A. cDFT calculations

To begin our cDFT calculations, first, we generate two bulk polymer solution phase diagrams for the two chain lengths we investigated, *N* = 10 and *N* = 40, as shown in Fig. 3. The phase diagrams allow us to determine the polymer–solvent interaction, *ɛ*_{ps}, which leads to either phase separation or complete miscibility across all polymer densities and guides the choices of parameters in the MD simulations. The phase diagrams are generated using analytic expressions for the free energy. Note that the phase diagrams in Fig. 3 are not typical thermodynamic phase diagrams because we have kept the values of *ɛ*_{pp} and *ɛ*_{ss} constant and only varied *ɛ*_{ps}; this is not the same as varying temperature, in which case all interaction parameters would change, but, instead, is just varying the solvent quality for the polymer. From the phase diagrams, we choose to use *ɛ*_{ps} = 0.93 kT and *ɛ*_{ps} = 1.0 kT as our poor and good solvent qualities for both polymer chain lengths. The pressure of the system is *p* = 0.650 909, which is the pressure of the pure solvent at a density high enough to remain completely liquid. All cDFT calculations are performed at this constant pressure as the bulk polymer density is varied by simultaneously varying the solvent density according to the relation $\rho sb=g(\rho pb;p)$ described in Sec. II; note that the relation between the two bulk densities depends on *ɛ*_{ps} and on chain length and so must be determined for each condition. Using the poor solvent, *ɛ*_{ps} = 0.93 kT, phase separation occurs at $\rho p,coexb=0.074\sigma \u22123$ and $\rho p,coex,2b=0.402\sigma \u22123$ for the *N* = 10 system and at $\rho p,coexb=1.02\xd710\u22124\sigma \u22123$ and $\rho p,coex,2b=0.565\sigma \u22123$ for the *N* = 40 system. Note that $\rho p,coexb$ is the polymer coexistence density on the polymer-poor side of the phase diagram, while $\rho p,coex,2b$ is the coexistence density on the polymer-rich side of the phase diagram.

Using the ALC algorithm implemented in Tramonto, we calculate the surface free energy and the excess adsorption of polymers at constant *ɛ*_{ps} as a function of bulk polymer density, starting at low polymer concentrations in the miscible region of the phase diagram. The surface free energy is defined as $\Omega s\rho \alpha r=\Omega \rho \alpha r\u2212\Omega b$, where Ω_{b} is the free energy of a bulk homogeneous reference system with the same density and composition. Figure 4(a) shows a set of Ω_{s} vs $\rho pb$ calculations performed for Δ*ɛ*_{wf} from 0.1 to 0.8. For Δ*ɛ*_{wf} > 0.1, Ω_{s} continuously decreases with increasing $\rho pb$ until the coexistence density, $\rho p,coexb$, is reached. Figure 4(b) shows the excess adsorption of polymers as a function of $\rho pb$ for all values of Δ*ɛ*_{wf}. The excess adsorption of species *α* is defined as

where *L* is the size of the cDFT computational domain and $\rho \alpha b$ is the bulk density of species *α*. As shown in Fig. 4(b), for Δ*ɛ*_{wf} > 0.1, the excess adsorption Γ_{p} diverges at the coexistence density, $\rho p,coexb=0.074$. This is the signature of a complete wetting transition that occurs as the bulk polymer density (or chemical potential) is increased at constant temperature toward $\rho coexb$. The monotonic increase in the divergence of Γ_{p} occurs for systems that are above the wetting temperature (in this case, above the wetting point determined by *ɛ*_{ps}/*kT*) and also above any prewetting transitions.^{69} These calculations thus suggest that at the bulk coexistence polymer density, a thick layer of polymers forms near the surface due to a complete wetting transition. This observation is consistent with previously reported wetting behavior in multi-component fluids and surface induced phase separation.^{69,70} Prior to complete wetting, increasing Δ*ɛ*_{wf} increases Γ_{p}, which means that more polymer segments are adsorbed onto the flat surface.

To confirm whether a complete wetting transition occurs with increasing polymer solution concentration at $\rho p,coexb=0.074$, the bulk phase separation density, we plot the density profiles of polymer segments as a function of distance from the surface. Figure 5 shows two plots of polymer density profiles for two bulk polymer densities before and after the wetting transition. In Figs. 5(a) and 5(b), polymer densities near the surface are perturbed from the bulk for more than 4 *R*_{g} (*R*_{g} = 1.3*σ* for *N* = 10) away from the surface for $\rho pb<\rho p,coexb$ (black solid line). Compared to Fig. 5(a) and 5(b) depicts the density profiles for a system with much higher Δ*ɛ*_{wf}. The amount of adsorbed polymers immediately adjacent to the surface is much greater, and the length scale over which the polymer density is different from the bulk density is also significantly longer when Δ*ɛ*_{wf} is higher. These results indicate that the surface interaction with the polymer influences polymer density over length scales (≈5*σ* for Δ*ɛ*_{wf} = 0.2 and ≈8*σ* for Δ*ɛ*_{wf} = 0.8) far beyond the range of the surface–polymer potential (*r*_{cut} = 2.5*σ*). Solvent density profiles corresponding to the two states in Fig. 5 are included in the supplementary material, Fig. S1. After the wetting transition, the polymer forms a (macroscopically) thick layer at the surface that extends to the entire computational domain. Note that the polymer density far from the surface in this case is equal to the coexistence density on the polymer-rich side of the phase diagram, $\rho p,coex,2b=0.403$. The *N* = 40 system also shows complete wetting transitions for Δ*ɛ*_{wf} > 0.4 but, interestingly, displays a first order wetting transition before reaching complete wetting for Δ*ɛ*_{wf} = 0.3; see the supplementary material for further details.

We next consider the case of a good solvent with *ɛ*_{ps} = 1.0. For sufficiently long chains, the adsorption of polymers onto an attractive flat surface from semi-dilute solution can be described by de Gennes’ self-similar “carpet” model.^{71–74} We note that, in general, our chains are too short to fit this model, although previous MD simulations did find the de Gennes scaling for sufficiently long chains.^{75} Figures 6(a) and 6(b) show Ω_{s} and Γ_{p} as functions of $\rho pb$ for *ɛ*_{ps} = 1.0. No wetting transition is observed as a divergence of Γ_{p} vs $\rho pb$ because the polymer solution does not phase separate at any bulk polymer concentration. When Δ*ɛ*_{wf} is small (Δ*ɛ*_{wf} < 0.3), Ω_{s} continuously increases as $\rho pb$ increases, while for Δ*ɛ*_{wf} ≥ 0.3, Ω_{s} decreases as $\rho pb$ increases with a gradual change in slope. The relationship between the surface excess adsorption of polymer, Γ_{p}, and $\rho pb$ also changes with Δ*ɛ*_{wf}. Note that when Γ_{p} < 0, there is less polymer at the surface than there would be if the polymer was at its bulk density throughout, while Γ_{p} > 0 indicates an excess of polymers over the bulk. We thus consider Γ_{p} = 0 to correspond to the adsorption (or desorption) transition. When Δ*ɛ*_{wf} < 0.3, Γ_{p} monotonically decreases as $\rho pb$ increases, but when Δ*ɛ*_{wf} ≥ 0.3, the relationship between Γ_{p} and $\rho pb$ is non-monotonic; Γ_{p} reaches a maximum value at small $\rho pb$ at a given Δ*ɛ*_{wf} and eventually becomes negative at a high $\rho pb$, where the adsorption–desorption transition occurs. The maximum in Γ_{p} implies the presence of an interplay between enthalpic and entropic driving forces. At lower bulk polymer densities (concentrations), the enthalpic attraction between the polymer and the surface drives the polymer to adsorb onto the surface. However, as the bulk polymer concentration increases, mixing entropy and the favorable interactions with the solvent reduce the polymer adsorption to the surface, until eventually the polymer desorbs.

To further examine the effect of surface–fluid interactions on excess polymer adsorption in the good solvent, we plot the maximum excess adsorption, Γ_{p,max}, as a function of Δ*ɛ*_{wf} in Fig. 7(a). For both chain lengths, increasing Δ*ɛ*_{wf} increases Γ_{p,max}, caused by more favorable surface–polymer interactions. The lower translational entropy of the long chains outweighs the configurational entropic penalty to confine the longer chains near the surface, so Γ_{p} is larger for the longer chains. To confirm that the increase in density near the surface is the cause of an increase in Γ_{p,max}, instead of a wetting transition, we plot the polymer density as a function of distance from the surface for all Δ*ɛ*_{wf} at the bulk polymer density where the maximum excess adsorption occurs for the *N* = 10 cases in Fig. 7(b). As Δ*ɛ*_{wf} increases, the locations of the peaks corresponding to monolayers packed near the surface and the far-field density remain constant, while the magnitudes of the density peaks near the surface are significantly increased. We quantify the amount of polymer segments adsorbed onto the surface as an areal density, Γ_{p,s}, by integrating the density profiles from the surface to the location of the first minimum in *ρ*_{p}(*x*) and converting to a 2D disk packing fraction,

where *l*_{min} is the distance from the surface to the first minimum in the density profile. In Fig. 7(c), we observe an increase in Γ_{p,s} as Δ*ɛ*_{wf} or N increases, similar to the trend in Fig. 7(a). The *N* = 40 system at Δ*ɛ*_{wf} ≥ 0.7 reaches and surpasses the random-close-packed limit for 2D hard disks, which is Γ_{p,s} ≈ 0.82–0.84. Thus, at the highest Δ*ɛ*_{wf} that is studied, the adsorbed polymer forms a dense layer that covers the surface. Together, the calculations in Fig. 7 show that the increase in Γ_{p,max} comes from the increase in surface coverage of the monolayers adjacent to the surface.

### B. MD simulations

Due to the large computational cost of cDFT calculations in 3D, we choose to use MD simulations to further our investigation of the effect of surface curvature on the polymer bound layer. We use the good solvent condition (*ɛ*_{ps} = 1.0) for the MD simulations to ensure that the polymers and the solvent do not phase separate, and we choose $\rho pb=0.1$, which is near the maximum adsorption from the cDFT calculations, and in the regime where we expect adsorption. After obtaining equilibrium configurations, we calculate the density profiles for each chain length, surface curvature, and Δ*ɛ*_{wf}, as shown in Fig. 8. The lack of layering peaks near the surface in the MD results is due to the roughness of the surfaces, which is neglected in the cDFT calculations. Systems with the smallest Δ*ɛ*_{wf} (black symbols) show a low density near the surface with density increasing monotonically from the surface to the bulk value, indicating that the polymers do not adsorb onto the surface. The overall trends in the density profiles show good agreement with the cDFT flat surface calculations with good solvent [Fig. 8(a)]. All perturbations in the polymer density calculated from MD simulations diminish within 6*σ* away from the surface, with slight increases as curvature and chain length increase. For Δ*ɛ*_{wf} ≥ 0.3, the polymer adsorbs onto the surface and the polymer density near the surface increases with increasing Δ*ɛ*_{wf} and with increasing N.

Previous theoretical and simulation work predicts that the polymer density profile near the flat surface should scale as *ρ*(*x*) ∼ *x*^{−4/3} in a central region between the adsorbed monolayer and the bulk density region.^{71–74} This central region should occur for distances larger than the cutoff in the polymer-wall attraction for *x* > 1.5*σ*. We show fits of the MD-calculated polymer density profiles for *N* = 40 to a power-law between 2 < *x* < 3.5 in Fig. 9. Our scaling region is very limited because of the short chain lengths in our simulations.^{75,76} For Δ*ɛ*_{wf} > 0.3, we observe good agreement between our simulation results and the theoretical prediction; our fits yield *ρ*_{p} ∼ *x*^{−1.32} for Δ*ɛ*_{wf} = 0.5 and *ρ*_{p} ∼ *x*^{−1.34} for Δ*ɛ*_{wf} = 0.7. Previous MD simulations of similar models found that *ρ*(*x*) ∼ *x*^{−4/3} for a somewhat larger region in *x* for *N* = 100.^{75}

Figure 10 shows the segmental and chain surface coverage as a function of Δ*ɛ*_{wf} for three curvatures and two chain lengths. Segmental surface coverage is defined as the number of polymer segments within 1.5*σ* of the surface per unit area (*σ*^{−2}), and chain surface coverage is defined as the number of independent chains with at least one monomer within 1.5*σ* of the surface per unit area. As Δ*ɛ*_{wf} increases, both the segmental and chain surface coverages increase in all cases. In Fig. 10(a), an increase in curvature or chain length both slightly increase the polymer segmental surface coverage, especially at high Δ*ɛ*_{wf}. In Fig. 10(b), the trend in the chain surface coverage as a function of Δ*ɛ*_{wf} is the same as that of the segmental surface coverage. Increasing curvature also increases the segmental surface coverage. However, unlike the segmental surface coverage, the chain surface coverage decreases with increasing chain length. In Fig. 10(c), the chain length difference is eliminated by multiplying the chain surface coverage by $Rg2=Nb/6$, where we assume that *b* = 1*σ* for simplicity, which measures the number of chains in a chain-size independent surface area. By doing so, the effect of chain length on the number of independent chains at the surface is reversed and matches that of the segmental surface coverage. Comparing Figs. 10(a) and 10(c), the normalized chain surface coverage shows a larger increase when the chain length increases than does the segmental surface coverage, meaning that the number of independent chains adsorbing onto the surface when the surface area is normalized by the polymer size responds more strongly to the increase in curvature. This shows that the segmental surface coverage is largely driven by the enthalpic interactions between the polymer segments and the surface, while the chain surface coverage is affected much more by chain lengths. The surface curvature also plays a role in the chain packing: longer chains tend to pack a higher percentage of their segments on highly curved surfaces, leaving fewer sites for another chain to occupy. The results in Fig. 10 suggest that while the enthalpic gain from forming contacts with the NP dominates the overall energetics, the translational entropy loss due to increases in chain length or surface curvature is non-negligible in determining the surface coverage of chains.

To further understand the energetics of segmental adsorption on the surface, we calculate the polymer segmental potential of mean force (PMF) from the density profiles, *w*(*r*)/*kT* = −log(*ρ*(*r*)), with *w*(*r*) = 0 as *r* → *∞*. We calculate PMFs using data from both cDFT calculations and MD simulations. Figure 11(a) shows the minimum values of *w*(*r*), *w*_{min}, as a function of Δ*ɛ*_{wf} for the flat surface and the two chain lengths for both cDFT and MD simulations. Although there is a slight difference in the values, the overall trends remain the same for both methods. When Δ*ɛ*_{wf} increases, *w*_{min} decreases, indicating that a larger amount of energy is required to remove a polymer segment from the surface. At the lowest value of Δ*ɛ*_{wf}, due to the small *ɛ*_{wp}, there is no energetic penalty associated with removing a polymer segment from the surface. However, as Δ*ɛ*_{wf} increases, the magnitude of *w*_{min} becomes far greater than that of *ɛ*_{wp}. For reference, the magnitude of polymer–particle interactions varies from *ɛ*_{wp} = 0.75 to *ɛ*_{wp} = 1.35. In addition, the effect of chain length on *w*_{min} depends on the Δ*ɛ*_{wf} value. At the lowest and at higher Δ*ɛ*_{wf} values, increasing chain length increases |*w*_{min}|, but at intermediate Δ*ɛ*_{wf} values (Δ*ɛ*_{wf} = 0.3), the chain length does not impact the segmental PMF potential well much for all curvatures. Figure 11(b) shows that for the 3D MD simulations with NPs, the trends are the same as in the flat surface case. For all combinations of Δ*ɛ*_{wf} and chain length, increasing curvature increases |*w*_{min}|.

Next, we investigate the adsorption dynamics of polymer segments on the surfaces. We identify segments that complete an entire adsorption process (adsorbing and desorbing) during the simulation, calculate the distribution of adsorbed time, and normalize the distribution so that the probability of adsorption at the smallest adsorbed time is *P*_{seg}(*t*_{min}) = 1. The resulting normalized distributions of segment adsorbed time are presented in Fig. 12. Figures 12(a) and 12(b) show the segment adsorbed time distribution for *N* = 10 and *N* = 40 systems, respectively. For clarity, only Δ*ɛ*_{wf} = 0.1 and Δ*ɛ*_{wf} = 0.7 calculations are included in Fig. 12; distributions for segment adsorbed time for all Δ*ɛ*_{wf} can be found in the supplementary material, Fig. S7. For both chain lengths and all curvatures, increasing Δ*ɛ*_{wf} increases the population of polymer segments that adsorb onto the surface for a longer time. This result is expected because the increase in polymer–surface interaction makes it more energetically unfavorable for the polymer to move away from the surface, which also agrees with the calculations of the PMF minima in Fig. 11. However, increasing the curvature of the surface, especially at Δ*ɛ*_{wf} = 0.7 (dashed lines, open symbols), leads to a minimal decrease in the segment adsorbed time. The small role of curvature is not surprising, given the short chain lengths studied; recall that *R*_{g}/*r*_{NP} only varies from 0.05 to 0.43 in our simulations, and longer polymers will be the subject of future studies.

Finally, we examine the polymer chain adsorption dynamics on the surfaces. Positions of independent chains are tracked during the MD simulation, and a chain is considered adsorbed to the surface if any of its segments is adsorbed (i.e., within 1.5*σ* of the surface). We use a similar procedure to normalize the distributions of chain adsorbed time so that *P*_{chain}(*t*_{min}) = 1. Figure 13(a) shows the chain adsorbed time distribution for all Δ*ɛ*_{wf} for the *N* = 10 flat surface systems. At Δ*ɛ*_{wf} = 0.1, the distribution appears to be Gaussian. However, as Δ*ɛ*_{wf} increases, not only do the distributions shift to longer times but also a large “shoulder” develops, suggesting two overlapping distributions of chain adsorbed times. Increasing Δ*ɛ*_{wf} increases the adsorbed time, which is consistent with prior research.^{77} Figure 13(b) shows the effect of curvatures on the chain adsorbed time distributions. At both Δ*ɛ*_{wf} shown in Fig. 13, increasing curvature shifts the distribution to shorter adsorbed times. The curvature effect is minor at low Δ*ɛ*_{wf}, but it becomes significant at large Δ*ɛ*_{wf}. The relatively small effect of the curvature for most of our systems is due to the small *R*_{g}/*r*_{NP}, which are not in the regime where the polymer size is comparable to the NP size. We anticipate that increases in *R*_{g}/*r*_{NP} could lead to larger differences in the adsorption dynamics. Although the effects due to curvature are relatively small compared to those from the polymer–surface interactions, it is evident that entropic effects that arise due to changes in packing resulting from changes in curvature are considerable in determining the polymer dynamics in the bound layer.

As mentioned in Sec. I, the reversibility of polymer adsorption is still a highly debated topic. Previous theoretical treatment of reversibility by Guiselin assumed that the chains in contact with the surface are irreversibly bound during solvent washing.^{78} Our simulations indicate that the adsorption of polymer chains is reversible and that the chains in the bound layer are actually quite mobile when exposed to a good solvent, which is supported by numerous experiments.^{79–83} The number of polymer segments adsorbed onto the surface remains relatively constant during the simulation time, but chains can desorb during this time, which suggests that chains can desorb a few segments at a time, while segments from a free chain can adsorb to occupy these previously occupied sites, consistent with previous studies.^{49,80,84,85} With larger Δ*ɛ*_{wf}, the chain adsorption times appear to consist of two overlapping populations. The tightly bound chains have a longer desorption time, and a separate population corresponding to the more loosely bound chains desorbs more quickly.^{47,84,86} This is consistent with the “parking problem” interpretation: the chains that adsorb first spread out to maximize the number of contacts with the surface, leaving little space on the surface for later chains to adsorb, which causes the “bimodal energy landscape” of chain adsorption.^{86,87} An additional complication is the aging of the bound layer that affects the desorption rate, which is not considered in this study.^{47,49} Without a doubt, more efforts are necessary to investigate these complications in the future.

## IV. SUMMARY

We use both cDFT and MD simulations to investigate the structure and dynamics of the bound polymer layer formed near a surface in solution. We vary the polymer–surface interactions, surface curvature, and polymer chain lengths. In the cDFT portion of this work, we place two polymer solutions with different solvent qualities next to a flat surface. In the poor solvent case, we calculate the density profiles below and above the bulk polymer concentration corresponding to the wetting transition and find that in the absence of a wetting transition, the density profile is a function of polymer–surface interactions and bulk polymer density. In the good solvent case, we find that the excess adsorption increases with increasing polymer–surface interaction strengths and polymer chain length. We perform MD simulations of polymer solutions in good solvent in contact with three surface curvatures (flat surface, large NP, and small NP). We examine the segmental and chain surface coverages and unsurprisingly find that an increase in surface–polymer interaction strength, surface curvature, or polymer chain length increases the amount of surface coverage by both polymer segments and chains. However, at first glance, systems with longer chains show lower chain surface coverage than systems with shorter chains, which is opposite to the trend in segmental surface adsorption. When normalized by the ideal chain radius of gyration, the longer chains show a larger number of independent chains adsorbing onto the surface, which is consistent with the segmental surface coverage. Finally, we examine the polymer segmental adsorption time between the bound layer and the bulk solution. The duration of a chain adsorption is most largely affected by the polymer–surface interaction strengths. The polymer chain adsorption time distributions suggest that while the polymer–surface interaction strength dominates the overall behavior, the surface curvature also impacts the dynamics at the chain level in a non-trivial way.

The following insights from our results will help guide experimental design. First, choosing the right solvent during the preparation of the PNC is of paramount importance to the thickness of the resulting polymer bound layer, echoing results found by Jouault *et al.*^{18} Second, the difference between surface–polymer and surface–solvent interactions is more important than chain length or surface curvature in determining the surface coverage of polymer segments and the dynamics of segmental adsorption, at least for short chains that are below the entanglement molecular weight and have *R*_{g}/*r*_{NP} < 0.5. The chain lengths we used in our study are not long enough to provide agreement with previous scaling arguments put forth by theoreticians, and future work involving longer chains is currently under way. Our study also suggests that in a polymer solution at higher temperature, adsorbed polymers are highly mobile and can desorb readily, which is consistent with the findings of Thees *et al.*^{47} Our molecular-scale investigation of the polymer bound layer shows that we can expect the adsorption of polymers to surfaces to be reversible under appropriate conditions although there may be other situations, such as rough “flat” surfaces, strong surface–polymer interactions, or higher *N*, where the adsorption is effectively irreversible.

## SUPPLEMENTARY MATERIAL

See the supplementary material for the details on the constant pressure density equations and solvent density profiles for *N* = 10; free energy, excess adsorption, and polymer density profiles for the *N* = 40 system using cDFT calculations; and the polymer segmental adsorption time distributions using the MD simulations.

## ACKNOWLEDGMENTS

The authors would like to acknowledge financial and computational support from the Department of Energy (DOE) under Award No. DE-SC0016421 and Sandia National Laboratories. This work was performed, in part, at the Center for Integrated Nanotechnologies, an Office of Science User Facility operated for the U.S. Department of Energy (DOE) Office of Science. Sandia National Laboratories is a multi-mission laboratory managed and operated by the National Technology and Engineering Solutions of Sandia, LLC, a wholly owned subsidiary of Honeywell International, Inc., for the U.S. DOE’s National Nuclear Security Administration under Contract No. DE-NA-0003525. The views expressed in this article do not necessarily represent the views of the U.S. DOE or the United States Government. We would also like to thank the DOE Science Graduate Student Research (SCGSR) fellowship for providing Emily Y. Lin an opportunity and financial support to work on this project in person at Sandia National Laboratories. In addition, the authors thank Dr. Eric Bailey for helpful discussions.

## DATA AVAILABILITY

The data that support the findings of this study are available from the corresponding authors upon reasonable request.