Interaction potential for coarse-grained models of bottlebrush polymers

Bottlebrush polymers are a class of branched macromolecules, consisting of numerous short chains densely grafted onto a linear polymer. They are rapidly emerging as attractive building blocks for soft materials that can be widely used in many applications.^1–9 Recent applications of this class of macromolecules include solvent-free soft elastomers,^10–16 surface modification for improved lubrication and antifouling,^17,18 photonic materials,^19–23 pressure sensitive adhesives,²⁴ self-healing polymeric networks,²⁵ and micellization for targeted drug delivery.^6,26 The wide adoption of the bottlebrush architecture stems from two distinct factors: (i) the presence of densely grafted side chains adds a non-trivial molecular thickness that reduces intermolecular entanglements and helps in rapid self-assembly^27,28 and (ii) many choices are available for the molecular design, e.g., degree of polymerization of the backbone and side chains,^29,30 grafting density,³¹ and monomer chemistry,³² which can be independently varied to tune material properties at the macroscopic level. Understanding the role of architectural parameters is crucial to improved designing of bottlebrush materials. In conjunction with experimental characterization and analytical theory, computational studies can provide key insights into the underlying physics³³ and help us correlate bulk material properties with the molecular architecture.

Bottlebrush polymer molecules are most commonly modeled using a coarse-grained bead–spring representation,^34–38 although several recent studies have used different models as well, such as density-functional theory,³⁹ self-consistent field theory,^40,41 and polymer reference interaction site models.^35,42 In these bead–spring models, each bead of the backbone and the side chains of bottlebrushes both typically resolve down to the level of a Kuhn step. However, resolving the side chains results in very large degrees of freedom at the molecular level significantly drives up the computational cost and makes it difficult to access the large length scale properties characteristic of multi-chain bottlebrush materials. For instance, a bottlebrush with 500 Kuhn steps on the backbone and 30 Kuhn steps/side chain with a grafting density of unity will require 15 500 beads for a single molecule. Clearly, this approach becomes restrictive for non-dilute bottlebrush solutions and melts, where the total number of molecules is large. This motivates us to seek a molecular model in which the side chains can be coarse-grained away in an appropriate fashion—we will call such models implicit side chain (ISC) models, in contrast to the more common explicit side chain (ESC) models. A typical simulation snapshot and schematic of the ESC model and a conceptual picture of the ISC model are shown in Fig. 1.

In our previous works,^43–45 we develop an ISC model for bottlebrushes based on a wormlike cylinder. This models a bottlebrush molecule using four descriptors: (i) a contour length L, (ii) a Kuhn length λ⁻¹, (iii) a diameter d, and (iv) an excluded volume parameter B.⁴³ These four values parameterize a bottlebrush viewed as a wormlike cylinder with a uniform diameter d rather than the individual architectural features, such as the length of the backbone or the side chains.⁴³ The parameters are determined from equilibrium averaged conformations obtained from molecular simulations using an ESC model and systematically fitted to analytical expressions for a wormlike cylinder.^43,46 The diameter d is approximately four times the gyration radius of the side chains, and the contour length L is found to be shorter than the actual contour length of the backbone due to local undulations within a tube populated by the side chain monomers.⁴³ For computational purposes, this model is discretized as a sequence of touching beads with an angle-dependent potential between successive bonded beads and a hard sphere like repulsive interaction between non-bonded beads.⁴³ This ISC representation was successful in capturing equilibrium and near-equilibrium properties, such as end-to-end distance, radius of gyration, hydrodynamic radius, and intrinsic viscosity, and accessing much longer bottlebrush lengths than is possible with an ESC model.⁴³ Furthermore, this same model is also able to accurately capture the stretching response of a bottlebrush under low force.⁴⁴ 

Rooted in these efforts, we applied this ISC model to study the self-assembly behavior of bottlebrush block copolymers in solutions.⁴⁵ We derived an intermolecular interaction potential from scaling arguments, with two parameters to be tuned to match experimental observations.⁴⁵ The intermolecular interaction is a soft potential, in contrast to commonly used hard-sphere potentials such as Weeks–Chandler–Anderson (WCA) potential.⁴⁷ Despite the fact that the simulation results based on the scaling-based potential are consistent with experimental measurements,⁴⁵ we would prefer to use the ESC model itself to determine the interaction between the coarse-grained beads of bottlebrushes using the ESC model. Meanwhile, with the derived potential in mind, we can directly compare the scaling relationships with simulations using the ESC model. It is beneficial for future studies to establish a systematic coarse-graining workflow so that similar ISC models for bottlebrush under different thermodynamic conditions (solvent qualities, monomer chemistries, etc.) can be generated.

In this paper, we investigate the interaction between two bottlebrushes and establish a protocol for incorporating these interactions into a coarse-grained, implicit side-chain representation of bottlebrush polymers. First, we use the ESC representation to calculate the potential of mean force (PMF), osmotic second virial coefficient, and interpenetration function between two bottlebrushes. We use a corresponding ISC representation that incorporates this molecular information (e.g., using the PMF to calculate interactions between coarse-grained beads, bending potential, etc.) and show that the two models give near-quantitatively matched PMF profiles and second virial coefficients. We note that we only consider the athermal condition in this work; however, our methods are not specific to this thermodynamic condition, and thus, the general methodology for understanding PMF between bottlebrush coarse-grained beads that we present here can be extended to different thermodynamic conditions. This can be useful for generating a coarse-grained potential of bottlebrushes, which is applicable in large-scale simulations of bottlebrush systems.

We consider two classes of models: (i) Explicit Side Chain (ESC) model, where each side chain is itself a bead–spring chain, and (ii) Implicit Side Chain (ISC) model, where the side chain degrees of freedom are coarse-grained out. Our goal is to calculate the PMF between bottlebrush polymer molecules using the ESC model and subsequently determine the interaction between the coarse-grained beads of the ISC model which generates the same intermolecular PMF, the latter being calculated based on an ISC representation alone.

Bottlebrush polymers molecules with explicit side chains are modeled by a set of beads connected by finitely extensible nonlinear elastic (FENE) springs. The backbone consists of N_bb beads where each backbone bead carries a side chain comprising N_sc beads [see Fig. 1(a)]. The motion of each bead i with position r_i is governed by a potential U for the system, which is given by

The potential consists of two contributions: the bonding potential U_s is between beads that are connected, and the excluded volume potential U_ev captures the pairwise repulsion between all beads. All beads interact via the purely repulsive WCA potential,⁴⁷ which depends on the distance r_ij = |r_i − r_j| between those beads,

where σ = 2a is the length parameter, ɛ = k_BT is the energy parameter (k_B is the Boltzmann constant and T is the absolute temperature), and r_cut = 2^1/6σ. In the case where bead i and bead j are connected, the form of the spring potential representing bonded interaction is

where k_s = 30ɛ/σ² and r_max = 1.5σ. This model has been used in our previous works^46,48 to study bottlebrush polymers with poly(norbornene) backbone and poly(lactic acid) side chains, where we demonstrated excellent agreement with experimental observations upon suitable parameterization for the chemistry-specific details.

The ISC model is a coarse-grained description of the ESC model that we discussed in Sec. II A. A detailed development of this model was presented in our earlier work,⁴³ and we adopt that overall formalism with several modifications to account for PMF data from the ESC model in the pairwise potential. Briefly, the ISC model treats a bottlebrush polymer molecule as a wormlike cylinder, characterized by a contour length L, Kuhn length λ⁻¹, diameter d, and an excluded volume parameter B.⁴³ These parameters are obtained by a systematic fitting procedure using conformational data determined from molecular simulations using the ESC model.⁴³ The wormlike cylinder is discretized into a sequence of beads of diameter d, connected by rigid bonds and arranged along the centerline of the molecule to generate a representation suitable for molecular simulations [see Fig. 1(b)]. The number of beads is proportional to the contour length of the backbone L_bb = l_b(N_bb − 1), where l_b = 1.94a is the equilibrium bond length in the ESC model.⁴³ By placing coarse-grained beads subsequently along a straight backbone, we have N_CGB = L_bb/d. However, the shrinking of the backbone due to local fluctuations further causes the contour length of the bottlebrush to be shorter than that of the backbone. This gives L = m_L(N_bb − 1), where m_L is a constant that accounts for this contour length disparity (also given in previous work⁴³). The bond length is thus l_CGB = L/N_CGB = m_Ld/l_b, and the fixed length of the bonds ensures that the contour length of an ISC chain does not fluctuate. By doing so, we obtain an integer number of coarse-grained beads to represent the corresponding ESC chain, while the subsequent beads are slightly overlapping. Thus, the total energy of the system using the ISC model is described by

To account for the molecular stiffness embodied in the Kuhn length, we use a bending potential of the form

where k_θ = (2λd)⁻¹ is the bending constant and $θij=cos−1ûj⋅ûi$ is the complementary bond angle between the unit vectors $ûi$ and $ûj$ along consecutive backbone bonds. We still need the pairwise interaction U_CGB between the coarse-grained beads to fully specify the ISC model. This will be extracted from the PMF data determined from the ESC model, as will be described in Sec. V.

To emphasize the key difference between the above described ISC model and that introduced in previous work,^43,45 our choice to represent the molecule with overlapping beads enables the incorporation of the PMF into the ISC model interaction potential (see Fig. 2 and discussion in Sec. II C). This is a consequence of how we set up the PMF calculation for the ESC model and later convert it to the interaction potential in the ISC model, noting that the shrinking of the bottlebrush contour is due to backbone fluctuations in the ESC model. We expect that this will help us better model bottlebrush polymers, especially in the case of inter-molecular interactions where the local shape of the molecule under the ISC representation is crucial.

We outline the workflow of this investigation in Fig. 2 to describe the interrelation between the various calculations we perform to both parameterize and demonstrate consistency between the ISC from the ESC models. We first use the ESC model to calculate PMF, shown in the left panel of Fig. 2 [including (a) and (b)]. This involves both calculating the molecular PMF in a traditional manner, where the molecule conformation is unconstrained [see Fig. 2(b)], and calculating the segmental PMF where the two bottlebrush molecules are aligned parallel to each other with straight backbones [see Fig. 2(a)]. In the latter case, the length of backbones are chosen carefully such that the ESC bottlebrush can be mapped to an integer number of coarse-grained beads, which is listed for various architectures in Table I. The segmental PMF can be normalized by this integer number to obtain the pairwise interaction between coarse-grained segments of the ISC model. Using this pairwise interaction as input, we then finalize the ISC model where the bond length is properly shortened, as displayed in Fig. 2(c). This enables us to sample the conformations purely based on the ISC model. We perform the molecular PMF calculation for the ISC model [see Fig. 2(d)] and relate it back to the initial point by comparing the results (PMF profiles, second virial coefficients, and interpenetration functions) with those from the ESC model. This comparison eventually serves as a metric to inform us the suitability of the ISC model that we build.

The PMF between two polymer molecules U(r) in dilute solution is defined as⁴⁹ 

where $U12r$ is the total intermolecular potential between two molecules with their centers of mass held at a distance r apart [see Fig. 3(a)] and the angle brackets denote an average over all possible orientations and all possible pairs of molecules. Computational studies using bead–spring models^50–54 typically evaluate $U12r$ by summing over the intermolecular interactions between the beads belonging to each molecule, i.e., without any intramolecular contributions; the operation indicated by the angle brackets is performed by averaging over multiple independent pairs of Monte Carlo (MC)-generated conformations and relative orientations. This method of direct evaluation via MC (which we call the single-chain method) is particularly suitable for athermal chains in dilute solutions, so we will use this technique for the molecular PMF between bottlebrushes as well.⁵⁵ For the segmental PMF calculation [see Fig. 3(b)], we adopt an alternative approach as a supplement—we perform dynamical simulations with two bottlebrushes present in the box and relax the side chain conformation simultaneously (which we call the two-chain method). In the long distance regime, the two methods provide identical PMF results (see Fig. S1 in the supplementary material); however, it can be difficult to sufficiently sample at short distances due to the strongly repulsive WCA potential and large overlap of the molecules. The benefit of two-chain method is that since the two molecules are relaxing at the same time, all these “heavily overlapping events” will be naturally avoided. This method also facilitates the extension of this workflow to more complex environments, such as polymer chains in the presence of ions or other components, as discussed by Withers et al.⁵⁶ 

In the molecular PMF calculation, we generated a set of equilibrium conformations using a combination of Brownian Dynamics⁵⁷ (BD) without hydrodynamic interaction augmented by global MC moves. The BD position update occurs following the discretized stochastic differential equation,

where R is the column vector of length 3N (N being the total number of beads) comprising all bead positions, U = U_s + U_ev is the total interaction potential, δt is the time step, and x is a column vector of length 3N containing random numbers drawn from a normal distribution that has a mean of 0 and a standard deviation of 1, which represents displacement due to Brownian motion. The right-hand side of the above equation was evaluated using the positions at time t and to obtain the bead positions at time t + δt on the left-hand side. In Eq. (7), the distance is non-dimensionalized with a, time with bead diffusion time τ = ζa²/k_BT (ζ is the bead friction coefficient), and energy with k_BT. Two different MC moves were used when necessary (see the supplementary material of Ref. 46): (i) backbone pivot and (ii) side chain double bridging, both subject to acceptance via the Metropolis criterion.⁵⁸ The BD time step was 10⁻⁴τ, with a pivot move every 50 BD steps and a side chain double bridging cycle every 11 BD steps. For every (N_bb, N_sc) pair, six trajectories were generated, each trajectory running for $>108$ BD time steps. This procedure generated $>6×104$ independent conformations on average for each bottlebrush molecule. An example of bottlebrush conformation using this model is shown in Fig. 1(a). The scheme we use to compute molecular PMF in Eq. (6) is as follows: We choose a pair of conformations from the set of equilibrium conformations and label them as molecule 1 and molecule 2. The center of mass of molecule 1 is then held fixed at the origin and that of molecule 2 is placed with its center at a distance r along the x-axis. Pairwise overlap between the beads of molecules 1 and 2 is checked using a hierarchy of axis-aligned bounding boxes.^59–61 For overlapping beads, we calculate the contribution of that pair toward $U12r$⁠. For each value of r, the relative orientation of molecule 2 is changed four different times by rigid body rotation to average out the angular dependence. However, this rotation is not necessary if the number of available pairs is large. After one sweep over all values of r for a given pair, the same procedure is performed over as many pairs as available. Note that in the molecular PMF calculation, no simulation box or periodic boundary conditions are used, as seen in Fig. 3(a). There is no dynamical evolution of the molecules and the calculation is purely based on the equilibrium conformations that are priorly generated. The above described strategy is also applied to both the ESC model and the ISC model for PMF calculation.

In the segmental PMF calculation between two aligned bottlebrushes, no MC moves are used, since the BD update for side chain beads is sufficient for generating independent conformations. The other aspects (e.g., BD time step and length of the simulation) are all kept the same as those used in the molecular PMF calculations. Two bottlebrush molecules were initialized inside a simulation box such that they cross the periodic boundaries along the y-direction, as seen in Fig. 3(b). In this case, we are essentially modeling two infinitely long and aligned bottlebrushes, calculating the interaction force between segments with a length dictated by the simulation box size in the y-direction. The backbone beads are fixed in a fully extended conformation, with the two backbones separated by a distance r in the x–z plane. Those two molecules are either using (i) two copies of configurations from the same configuration set or (ii) two independent molecules that are placed in the box simultaneously; these correspond to PMF calculation methods 1 and 2, as previously described. We vary N_sc from 2 to 20 and choose N_bb to be such that the corresponding number of coarse-grained beads N_CGB is close to an integer number. The chosen values of N_bb as well as the corresponding N_CGB can be found in Table I. The box size along y is set to be L_y = N_bbσ, which ensures that the backbone bond crossing the periodic boundaries along this direction has the same length as all other bonds, satisfying the periodic boundary condition. The box size along x and z directions is taken to be large enough to prevent interaction between the periodic images in those directions. Typically, the box extends by more than ten times the side chain radius of gyration beyond the two molecules along x and similarly along z directions. The average force between two bottlebrush molecules at each distance is obtained after a sufficiently long simulation, and the segmental PMF is calculated via numerically integrating the force using trapezoidal rule. The discretization of the distance is small enough to ensure a good numerical approximation from this integration. Beyond the longest distance that is calculated, the force will be linear extrapolated to obtain a residual triangular area, which contributes a small amount of energy to the overall PMF profile.

We focus exclusively on a segmental PMF for aligned backbones to maximize the points of interaction between bottlebrush molecules to benefit the sampling of the PMF. It is also possible to consider unaligned bottlebrush backbones, for example, with the two backbones positioned orthogonally. For computational simplicity, we do not perform this calculation; in part, we justify this by noting that the difference between a parallel vs orthogonal orientation is still represented even in the ISC due to the connectivity between beads. We also argue based on the results in Sec. V that this is ultimately unnecessary, as pair potentials calculated from aligned bottlebrush PMFs give rise to excellent matching between ISC and ESC models.

Once molecular PMF U(r) is obtained, the osmotic second virial coefficient A₂ was obtained using the formula⁴⁹ 

where N_A is the Avogadro number, ⟨S²⟩ is the mean square radius of gyration, M is the molecular weight, and $ρ=r/⟨S2⟩1/2$⁠. For the bottlebrush molecules studied here, M = Nm₀, where m₀ is the molecular weight of a segment and $N=Nsc+1Nbb$ is the total number of beads. Note that ⟨S²⟩ can be directly determined from the set of equilibrium conformations. The interpenetration function ψ was obtained from A₂ using the formula⁴⁹ 

We first focus on the results for the PMF between unconstrained bottlebrush polymers. The benchmark results of linear polymers are shown in Fig. S2 of the supplementary material. We highlight two example series of data in Figs. 4 and 5. First, we fix the backbone length as N_bb = 320 and vary the side chain length from N_sc = 0–32, as shown in Fig. 4. Here, the case with N_sc = 0 serves as a linear polymer baseline for comparison. As the side chain length increases, the PMF profile exhibits stronger repulsion between bottlebrush molecules. The upward shift in the PMF is mainly induced by the increased steric repulsion between side chains from different bottlebrushes, which is a consequence of increasing side chain length. However, the extent of the PMF is well-described at all side-chain lengths N_sc by the radius of gyration $⟨S2⟩1/2$⁠, suggested by the converging behavior of the curves at $r/⟨S2⟩1/2≈2.0$⁠.

Second, we fix the side chain length as N_sc = 14 and vary the backbone length from N_bb = 20–640, shown in Fig. 5(a). In this case, the bottlebrush conformation is closer to a spherical object at the low backbone length limit, while at the long backbone length limit, it acts more like a linear polymer. Therefore, when normalized by $⟨S2⟩1/2$⁠, the length scale is quite different among different N_bb. With N_bb = 20 and N_bb = 40, at the distance $r/⟨S2⟩1/2=0.8$⁠, the two bottlebrush molecules overlap significantly because the side-chains must interpenetrate due to their star-like structure. This gives rise to a marked increase in the PMF due to the lack of acceptable configurations. To contrast, for N_bb > 80, the PMF profile resembles those in Fig. 4 and goes to a finite value of U(r) because the coiled structure of the bottlebrush at these lengths allows them to overlap without significant side-chain interpenetration. To evaluate this change in molecular shape in a quantitative manner, we use the volume fraction N/V within pervaded volume $V∼⟨S2⟩3/2$ as a metric relevant for bottlebrush–bottlebrush overlap. As can be seen in Fig. 5(b), N/V generally decreases as the backbone length increases, indicating a transition of the bottlebrush conformation from a spherical object to a linear-like chain. This is also concomitant with the change in the PMF profile, where the shape of the PMF profile resembles that of the linear polymer in the long-chain limit, while having a hard sphere-like shape in the short-chain limit.

Figure 6 shows the second virial coefficient A₂ for linear and bottlebrush polymers as a function of backbone and side chain lengths. For linear chains, scaling theory predicts⁶² A₂ ∼ N^3ν−2, which for athermal chains (ν = 0.588) reduces to A₂ ∼ N^−0.236. Fitting to our data for A₂ for long chains, we obtain a scaling exponent of −0.29. This value is slightly lower, indicating that our simulations are still not in the long-chain limit, even for N_bb = 1280. To compare, earlier work by Dautenhahn and Hall⁵² reported an exponent of −0.272 and Yethiraj and Hall⁵⁰ reported −0.25 after extrapolating their data to N → ∞. Our prediction is in good agreement to these earlier results. Compared to linear chains, the second virial coefficient of bottlebrushes are systematically lower, with A₂ decreasing with an increase in N_sc with N_bb held fixed. We attribute this to the increased size of the bottlebrush chain as the side-chain interactions induce an increasing molecular rigidity, leading to less compact structures. This is also consistent with the general trend of decreasing A₂ with an increase in N_bb in all curves.

Figure 7 shows the interpenetration function ψ as a function of overall and backbone molecular weight for both linear chains and bottlebrush polymers. This quantity, as suggested by the name, is a good measurement of how much two molecules are interpenetrating each other at a similar normalized distance (with respect to radius of gyration). For sufficiently long molecules, ψ of both linear chains and bottlebrush polymers can be seen to approach an asymptote. Renormalization group theory⁶³ predicts that for a self-avoiding chain, ψ → 0.219 as N → ∞. We show this as a dashed line in Fig. 7, noting that our results are indeed consistent with this asymptoting value. In the case of bottlebrushes, there is a non-monotonic trend in ψ with a peak at short backbone lengths—this reflects the change in the overall shape of a bottlebrush polymer molecule from a spherical star to a cylindrical rodlike shape before transitioning to a chain-like conformation as the backbone length is increased at a fixed side chain length. This trend is related to our observations connecting the intrachain monomer volume fraction to the molecular shape, in that the interpenetration is a way to quantify the effect of short-range interactions on bottlebrush–bottlebrush interactions. The results of segmental PMF between aligned bottlebrushes are shown in Fig. 8, with both the original data from BD simulations (two-chain method) and a series of fitting based on scaling arguments. For all the side chain lengths we examine here, two-chain method provides a rather smooth curve throughout a large range of length scale. We show in Fig. S1 of the supplementary material that the single-chain method is quantitatively matching but becomes unreliable at the low distance regime due to the inability to adequately sample highly interpenetrating conformations. To better understand the trend in these segmental PMF profiles, we apply the scaling analysis that we developed in our previous work to fit the data to an analytical form.⁴⁵ Here, we adopt the following expression:

where the exponent $m=13ln⁡28ln(2/d)$⁠. The constant C includes all molecular information that is not a function of N_sc, as well as the proportionality constant inherent to the scaling analysis that was used to develop this expression. We choose d to be the point at which the PMF goes to $∼0.1kBT$ (0.2k_BT for N_sc = 20) to reflect that the scaling argument does not account for long-range fluctuations in interactions due to side-chains that can “stretch” beyond d. We thus expect long interaction tails, which we observe as the difference between the scaling theory and simulation data at large N_sc. We choose C = 3.2 as an overall fitting parameter of order unity that leads to excellent matching over all side-chain values N_sc (Fig. 8). For short side chains (N_sc = 2 and N_sc = 4), there is a small degree of deviation that we attribute to the inherent limitations of the scaling picture in describing short chains.

The results we obtain here suggest that the scaling-based theory we derived in our previous work⁴⁵ is already a reasonable approximation of inter-bottlebrush interactions, and in particular, the simple form of this potential is enabling for large-scale simulations that mainly focus on materials properties at a large length scale, such as block bottlebrush self-assembly. This provides a connection between this current study of segmental PMFs from the ESC model and experimentally observable properties; the results in this paper are consistent with pairwise interactions that are used to predict scattering and phase behavior in block bottlebrush polymers.⁴⁵ Now, we have gone further and can show how to incorporate ESC-determined PMFs directly into the ISC model without needing to assume a theoretical form to make predictions of bottlebrush solution properties.

We can use the PMF as an effective pairwise potential in a coarse-grained ISC representation of a bottlebrush so that it is informed by the ESC calculations. We fit a power-law curve to the segmental PMF data from Sec. IV and extrapolate this to extend over all values of r. This PMF can then be rescaled and directly incorporated as a pair potential into the ISC model, for which standard MC is performed to generate equilibrium bottlebrush conformations.

The PMF from both the ISC and ESC model can be compared by establishing the correspondence between the backbone length and the number of coarse-grained beads. The PMFs for both ISC and ESC models are plotted in Fig. 9 for N_sc = 14 using the N_CGB values listed in Table I to find the values of N_bb associated with the ISC representations. We include analogous plots of N_sc (Fig. S3) in the supplementary material. In general, for all N_sc we considered, we observe near-quantitative matching between the ISC and ESC models. Differences are most noticeable when the backbone is relatively short (i.e., small N_CGB) in the ISC model. This is attributed to two reasons: (i) end effects, which were not part of the ESC segmental PMF parameterization, become significant in this limit, and (ii) the value of $⟨S2⟩1/2$ that we use as a normalization constant for the x-axis is least accurate in this same limit. This latter point is rationalized by the use of only the center of the coarse-grained beads as a proxy for the mass distribution in these chains, despite the non-negligible distribution of mass throughout the coarse-grained bead. While these deviations are noticeable, they mostly become negligible at longer bottlebrush lengths where there is a clear consistency between these two approaches. Indeed, the main motivation to develop the ISC model is to simplify the calculation of high molecular weight bottlebrushes, which is where the agreement between ISC and ESC PMFs is the most quantitative.

We can further consider a comparison of the second virial coefficients and interpenetration functions between ISC and ESC models. As shown in Fig. 10, A₂ from the ISC model is systematically smaller than that of the ESC model. This is mainly attributed to the error in the calculation of the molecular size, which will be enhanced when $⟨S2⟩3/2$ enters the calculation. This disparity is largely normalized out in the PMF profiles such that the integral in Eq. (8) is essentially identical between the two models, and the disparities primarily arise in the prefactor. Thus, the systematic offset in A₂ indicates that one of the biggest disadvantages of the ISC model is that the size of the molecule is underestimated, but this would be a general consequence for all other similar coarse-grained models. In addition, the ISC calculation of A₂ is not monotonically decreasing at the short-chain limit, which we attribute to the aforementioned inaccuracies at the short-chain limit. This is even more notable when calculating the interpenetration functions, plotted in Fig. 11, where ψ shows a significantly larger peak at small values of N_bb. Notably, however, these significant differences become quite small at the long-chain limit where this approach works well.

In this paper, we studied the PMF between bottlebrushes using both an ISC and the ESC model. The molecular PMF profile between bottlebrushes shows a similar trend as the profile of a linear polymer when the bottlebrush is long but shows a qualitative difference when the bottlebrush is short. This trend in the PMF profile is, in part, related to the shape change of a bottlebrush molecule, where it behaves like a spherical object when the backbone is short but resembles a linear chain when the backbone is long. We calculated the segmental PMF between aligned bottlebrushes and found that this PMF can be described by the previously derived scaling form in a precise manner. Finally, by mapping the segmental PMF between aligned bottlebrushes to a potential between coarse-grained segments, we incorporated ESC information directly into an ISC bottlebrush model. We used this ISC model to perform a molecular PMF calculation and found near-quantitative agreement between the ISC and the ESC models. We note that the ISC model is most inaccurate when bottlebrushes have short backbones; however, this disadvantage occurs in the portion of the parameter space where the advantage of the ISC model—the reduction of computational cost vs the ESC model—is not as significant. There are a number of ways in which this model could be improved to address this limitation for short bottlebrushes, for example, explicitly parameterizing the chain ends interactions in the ESC model, although this comes at additional computational cost. For longer bottlebrushes, where the ESC model is particularly costly, the incorporation of these PMF calculations is a useful way to inform much more efficient ISC models that are useful for simulating multi-chain bottlebrush materials. In the future, this methodology holds promise as a way to include other useful interactions (charges, crosslinks, χ-driven interactions, etc.) informed by chemistry considerations.

See the supplementary material for (i) a comparison of the single-chain and two-chain methods for calculating the PMF, (ii) PMF calculations for polymers with different lengths, and (iii) comparisons of the ISC and ESC PMFs for a wide variety of side chain and backbone lengths.

This work was supported by the National Science Foundation under DMREF Award No. DMR-1727605 and the ACS Petroleum Research Fund under Award No. 61500-ND7.

The authors have no conflicts to disclose.

T.P. and S.D. contributed equally to this work.

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