Simulation of weak polyelectrolyte brushes: the effects of ionizable monomer fraction and monovalent salt

Abstract

We perform simulations using the grand reaction method (GRM) to investigate the effects of pH, monovalent salt concentration, and ionizable monomer fraction f on the swelling of weak polyelectrolyte brushes (PEBs) consisting of polymer chains with random sequences of neutral and weakly basic monomers. The brush height responses of the model PEBs obtained by incrementally increasing or decreasing pH are found to be indistinguishable. This behavior sharply contrasts with the pronounced hysteresis observed in some experiments, indicating that these measurements probe long-lived metastable states. Upon increasing salt concentration, the model PEBs exhibit a crossover from the osmotic brush regime to the salted brush regime, as predicted by mean-field theory. However, as observed in experiments, the scaling exponents in both regimes depend sensitively on f and pH and are found to be generally smaller in magnitude than predicted. Nonetheless, we demonstrate that the expected scaling behavior can be approximately recovered under specific conditions. These results support the growing body of work elucidating the complex behaviors of stimuli-responsive polymer brushes and highlight outstanding challenges in reconciling findings from experimental, simulation, and theoretical studies.

---

1 Introduction

Polyelectrolyte brushes (PEBs) are composed of charged polymer chains that are anchored by one end to a solid surface. They modify the interfacial properties of materials, imparting new stimuli-responsive functionality that can be tailored for applications ranging from actuation¹ and sensing^2–5 to drug delivery⁶ and separations.^7–9 The interfacial properties of PEBs are largely regulated by their conformations, which are in turn influenced by the polymer molecular weight, grafting density, and charge distribution. Controlling the interfacial response thus requires understanding how these key parameters affect PEB conformations under different conditions. Nonetheless, the mechanisms underpinning the responses of PEBs to changes in critical environmental factors such as solution pH and ionic strength remain incompletely understood.

Depending on their response to environmental conditions, PEBs may be classified as either weak or strong. The protonation state of weak PEBs varies reversibly in response to changes in pH and salt concentration, whereas strong PEBs remain fully ionized across a wide range of solution conditions. Mean-field theories predict several distinct conformational regimes for weak and strong PEBs.^10,11 Specifically, they predict that PEBs transition from the osmotic regime to the salted brush regime as the concentration of monovalent salt is increased above the crossover concentration . These regimes are distinguished by the scaling behavior of the brush height h. In the osmotic regime, the brush height for weak PEBs is predicted to increase with c_salt, scaling as , where α_b is the degree of ionization in the bulk solution, N_m is the number of monomer segments, and ρ is the grafting density. For strong PEBs, h ∼ f^1/2N_m increases with the fraction of charged monomers f but is independent of c_salt. In the salted regime, by contrast, h is predicted to decrease with c_saltvia h ∼ β^2/3N_mc_salt^−1/3ρ^1/3, where β is α_b and f for weak and strong PEBs, respectively.

These mean-field theory predictions have been tested experimentally. For strong PEBs, experimental measurements follow theoretical scaling laws, where brush height depends on N_m and f in the osmotic regime and N_m, f, ρ and c_salt in the salted regime.^12–16 For weak PEBs, tests of the scaling predictions have proved more challenging. The brush height is predicted to vary non-monotonically with increasing monovalent salt concentration as the system crosses over from the osmotic to salted brush regimes. This predicted behavior is qualitatively consistent with that observed in experiments. However, analysis of the dependence of h on salt concentration yields scaling exponents that vary strongly with solution conditions and brush structural parameters and generally have smaller absolute values than predicted in both the osmotic and salted brush regimes.^17–25 Furthermore, this non-monotonic behavior is also observed under experimental conditions where the weak PEBs are highly charged and expected to behave like strong systems. Thus, there remains ambiguity regarding the conditions where strong versus weak behavior is expected.

Molecular simulation studies have also been performed to understand the behavior of PEBs. Strong polyelectrolytes can be readily simulated using standard fixed charge models. Simulation studies of strong PEBs are largely consistent with mean-field theory predictions and experiments.^26,27 Simulating weak polyelectrolytes, however, requires employing advanced computational methods to model charge regulation and predict changes in the polymer's ionization state with varying environmental factors, such as pH and salt concentration.^28–33 As a result, comparatively few molecular simulation studies of weak PEBs have been performed.^34–38 Thus, simulation studies on weak PEBs offer new opportunity to scrutinize scaling theory predictions and understand the molecular origins of their apparent discrepancies with experiments.

Here, we use molecular simulation to investigate the responses of coarse-grained models of weak PEBs to changes in pH and monovalent salt concentration. Inspired by recent experiments, we examine PEBs consisting of polymer chains with random sequences of neutral and weakly basic monomers and vary the monomer composition to study its effects on brush response. Charge regulation in the PEBs is simulated using the grand reaction method (GRM). The GRM rigorously accounts for the effects of solution conditions and local environment on monomer ionization, and it has been applied to model charge regulation in solutions,^39,40 brushes,^35–38,41 and networks^42–44 of ionizable polymers. The GRM simulations successfully predict a crossover between osmotic and salted brush regimes. We characterize the scaling behavior of the brush height with monovalent salt concentration in each regime, as well as the dependence of the crossover salt concentration on polymer composition. Consistent with experiments,^24,45 the scaling exponents from the GRM simulations in both regimes are found to be generally smaller in magnitude than predicted by mean-field theory. Nonetheless, we are able to identify conditions where reasonable agreement with theory is obtained. Although our findings help to reconcile observations from experiments, theory, and molecular simulation, they also highlight several outstanding issues that present challenges in making direct comparisons between them.

2 Methods

2.1 Models

We performed simulations using ESPResSo 4.2⁴⁶ to investigate the effects of monovalent salt concentration and the ionizable monomer fraction on the pH response of planar weak PEBs (Fig. 1). The models and method are similar to those used in our recent studies of weak polyampholyte brushes.^41,47 Polymers were modeled as coarse-grained, bead-spring chains with N_m = 80 monomers each. Monomers were treated as single beads and assigned an identity of either a weakly basic (B) or neutral (N) chemical moiety. Titration reactions involving the weakly basic monomers (B ⇌ HB⁺ + OH⁻) were simulated using the grand-reaction method (GRM).^42,43,48,49 The solvent (water) was treated implicitly as a dielectric continuum, and the salt species (Na⁺ and Cl⁻) and other free ions (H⁺ and OH⁻) were modeled as charged beads. Model parameters and physical quantities from the simulations are reported using m, σ, and k_BT as the fundamental units for mass, length, and energy, respectively, where k_B is the Boltzmann constant and T is temperature.

Fig. 1 Model PEB with ionizable monomer fraction f = 0.6. The dimensions of the system are L_x = L_y = 46.2σ and L_z = 160σ. Grey surfaces are the lower and upper walls at z = 0 and z = L_z, respectively. The brush is in equilibrium with a bulk reservoir with pH^res = 7 and salt concentration c^res_salt = 0.1 M. The dearth of H⁺ and OH⁻ species is due to the high salt concentration of the system.

Planar brush configurations were created by end-grafting M = 64 polymer chains to a surface in the x–y plane of the simulation cell by fixing the position of one of the free ends at z = 1σ on a square lattice with areal density ρ and spacing d = ρ^−1/2 (Fig. 1). The x × y × z dimensions of the simulation cell were L_x × L_y × L_z, where L_x = L_y = L = (M/ρ)^1/2 and L_z = 2N_mσ. The PEBs were simulated using a slab geometry that was created by placing bounding lower and upper walls at z = 0 and z = L_z, respectively, and imposing periodic boundary conditions along the x and y directions. The chemical composition of the polymers was specified by generating random sequences of B and N type monomers to produces chains with the targeted ionizable monomer fraction, f = N_B/N_m, where N_B is the number of basic monomers.

Bonds between neighboring monomers along each chain were modeled by the finite extensible nonlinear elastic (FENE) potential,^50,51

where r_ij is the scalar distance between particles i and j, Δr_max = 1.5σ is the maximum bond length, and k_bond = 30k_BT/σ² is the bond stiffness.

Excluded volume interactions between all particles (polymer monomers and free solution ions) were described using the Weeks–Chandler–Andersen (WCA) potential,⁵²

Length scale and interaction strength parameters of σ_ij = σ and ε = k_BT, respectively, were used to give all particles identical excluded volume interaction.

Lastly, a Coulomb potential was used to model electrostatic interactions between charged species,

where z_i and z_j are the charge numbers of particles i and j. A Bjerrum length of λ_B = 2σ was used in the simulations. This choice implies a characteristic length scale of σ ≈ 0.355 nm because λ_B ≈ 0.71 nm for water at room temperature. Long-range contributions to the electrostatic interactions were treated using a three dimensional particle–particle–particle–mesh (P3M) solver. The electrostatic layer correction (ELC) method was used with a standard gap size of L_ELC = 0.2L_z to account for the two-dimensional slab geometry of the PEBs.^53–55 Other parameters, including the cutoff for the Coulomb potential, were chosen to ensure a relative error of 10⁻⁵ in the computed electrostatic forces.^53–55

2.2 Acid–base equilibrium

The GRM was used to model acid–base equilibrium in PEB systems coupled to an external reservoir with free solution ion concentrations c^res_i for species i ∈ {Na⁺, Cl⁻, H⁺, OH⁻}. The following reactions were simulated for the weakly basic species (B):where pK^base_A is the intrinsic dissociation constant of the base, K_w = 10^−pK_w is the ion product for water, and the other Ks denote solution-phase reaction equilibrium constants.

Additionally, the following auxiliary reactions were performed to maintain equilibrium with the reservoir:

where is the ionic strength, μ_ex(I^res) is the excess chemical potential for a neutral ion pair, and ∅ is an empty set associated with insertion/deletion of the neutral ion pair.

As discussed in ref. 41 and 48, the reaction network (eqn (4)–(10)) couples the system to an external bulk reservoir and allows simulation of acid–base equilibrium over a broad range of solution conditions. In our simulations, the physical parameters pK_w = 14 and T = 1/k_B are fixed, and we vary the reservoir conditions (pH^res = 3–11 and c^res_salt = 10⁻⁶–10⁰ M) and the ionizable monomer fraction f = 0–1 to study their effects on PEB behavior. An intrinsic dissociation constant pK^base_A = 8.4 and grafting density ρ = 0.03σ⁻² were chosen to mimic the systems examined in recent experimental studies in which ionizable 2-(dimethylamino)ethyl acrylate (DMAEA) and neutral 2-hydroxyethyl acrylate (HEA) monomers were used to create PEBs with different f.^45,56 The excess chemical potential μ_ex(I^res) in eqn (10) was calculated from molecular dynamics (MD) simulations of the bulk reservoir using the Widom test particle method^57,58 and procedures identical to those reported in ref. 41 and 48.

2.3 Sampling protocol

The properties of the model PEBs were sampled using a combination of molecular dynamics (MD) and Monte Carlo (MC) moves, as described in detail in our previous study.⁴¹ Briefly, we employed two types of MC moves. In a reaction MC move,^48,59,60 one of the forward or reverse reactions in eqn (4)–(10) is selected at random and then subsequently attempted by replacing randomly selected reactant particles with the corresponding products. In a charge exchange MC move, the identity of {B, HB⁺} monomer pairs are swapped while keeping the number of each type of species and total charge constant.^61–63 The latter type of MC move expedites equilibration of the charge distribution along the polymer chains.^61–63 The Metropolis-Hastings acceptance criteria for these MC moves are reported in ref. 41.

Sampling was conducted by performing MC-MD cycles, with each cycle consisting of one MC and one MD step. Each MC step consisted of 2N_m attempted reaction moves and N_m/10 attempted charge exchange moves. Each MD step consisted of 1000 time steps using the velocity Verlet integration algorithm (time step δt = 0.01τ, where τ = σ(m/k_BT)^1/2) and a Langevin thermostat (friction coefficient γ = τ⁻¹). Equilibration of the PEB systems was performed for 10⁴ MC-MD cycles, followed by a production phase of equal duration. Time series (autocorrelation function) analysis showed that the production phase of each simulation was sufficient to generate more than 50 statistically independent samples for each observable. Statistical uncertainties for each observable were estimated from the standard error of the mean. In each case, the reported results were obtained from simulations of a single PEB system containing M = 64 chains with unique random sequences but the same ionizable fraction f. Finite size effects were assessed by comparing the results with larger systems and multiple independent replicas with different realizations of the random chain sequences. Agreement between the datasets demonstrated that a single replica with M = 64 chains was sufficient to represent the ensemble of random chain sequences with a given f and that the simulations were adequately sampled using the procedures described above.

3 Results and discussion

We first characterize the ionization state of the PEBs as a function of pH^res (Fig. 2). The degree of ionization is computed as the average fraction of basic charged monomers:

Fig. 2 Titration curves (average degree of ionization 〈α〉 vs. pH^res) for brushes with different ionizable monomer fractions f and salt concentrations of c^res_salt = (a) 1 × 10⁻⁶ M, (b) 1 × 10⁻² M, and (c) 1 × 10⁻¹ M. Ideal Henderson–Hasselbalch titration curves (eqn (12)) are shown as dashed lines. Uncertainties are smaller than the symbol size.

In the ideal limit, the variation of this quantity with pH is given by the Henderson–Hasselbalch (HH) equation:

where pOH = pK_w − pH, and pK_B = pK_w − pK^base_A. The titration curves for the PEBs (〈α〉 vs. pH^res) exhibit sigmoidal shapes and decrease monotonically as pH^res increases (Fig. 2). Although this behavior is qualitatively consistent with that predicted by the HH equation, significant quantitative deviations are observed. These deviations result from the Donnan partitioning of free solution ions and non-ideal behavior arising from the presence of strong electrostatic repulsions between the ionized monomers.⁴⁸ As the repulsions are enhanced, ionization becomes less favorable at a given pH^res, resulting in larger deviations from ideality. Accordingly, the deviations from the HH equation become increasingly pronounced as the ionizable monomer fraction f increases. The deviations also increase as the salt concentration c^res_salt decreases because of a reduction in the screening of the electrostatic repulsions along the polymer chains by free solution ions. This observation is consistent with experiments on weak basic PEBs, which show that ionization is hindered upon decreasing salt concentration, as reflected by a decrease in the apparent pK_a extracted from the titration curves.^45,56,64

To examine changes in the structure of the PEBs, we calculate brush height,²⁶

where ρ_m(z) is the density of monomers as a function of the distance from the grafting surface z (Fig. S1 in SI). The factor of 2 in eqn (13) is often omitted in the definition of brush height but has been included here to ensure that 〈h〉 approaches the contour length in the limit that the polymer chains become fully stretched.

Similar to the trend observed in the titration curves, the brush height exhibits sigmoidal-like behavior and decreases monotonically with increasing pH^res (Fig. 3). At high pH^res, where 〈α〉 →0, 〈h〉 for the ionizable PEBs (f > 0) saturates at values similar to that of the neutral brush (i.e., f = 0). At low pH^res, by contrast, 〈α〉 → 1 and the electrostatic repulsion along the chains cause the ionizable PEBs to swell. The degree of swelling at low pH^res increases with f due to the enhanced electrostatic repulsions. Additionally, the onset of swelling shifts to higher pH^res as c^res_salt increases due to increased ionization under these conditions (Fig. 2; Fig. S2 in SI).

Fig. 3 Brush height 〈h〉 vs. pH^res for PEBs with different ionizable monomer fractions f and salt concentrations of c^res_salt = (a) 1 × 10⁻⁶ M, (b) 1 × 10⁻² M, and (c) 1 × 10⁻¹ M. Uncertainties are smaller than the symbol size.

Interestingly, experimental swelling measurements on PEBs often exhibit pronounced hysteresis and history-dependent responses to changes in pH^res and relative humidity, indicating that the brushes become trapped in long-lived metastable states.^65–67 To check for this behavior, we performed simulations in which pH^res is cycled through incremental step changes. In sharp contrast with experiments, we observe no differences between the pH sweeps performed in the increasing (3 → 11) or decreasing (11 → 3) directions (Fig. 4), indicating that our simulations sample reversible equilibrium behavior. Charge hysteresis has been observed for PEs in solution, where it has been attributed to charge-conformation coupling in which compact chain conformations resist ionization more strongly than extended ones.^68,69 Similar charge-conformation coupling mechanisms have also been invoked to explain the swelling hysteresis in PEBs.^65,70,71 It has been proposed, for example, that the swelling hysteresis in PEBs may arise from the free chain ends near the periphery being more easily ionized due to their greater solvent accessibility.^45,56,70 Additionally, it has been posited that monomer hydrophobicity causes swelling hysteresis by hindering transport of bulk solution into the brush and inhibiting ionization of deeper brush segments.^{23,56,65,72,73} The hypothesized importance of water, hydrophobic effects, and transport into the brush suggests that there are several possible reasons that the simulations do not capture the experimentally observed hysteresis behavior. The implicit treatment of the solvent as a uniform dielectric continuum prevents the PEB model from describing water's structural properties (e.g., hydrogen bond networks) and its inhomogeneous distribution across the interface and within the brush phase. The PEBs are modeled under good solvent conditions, which are mimicked by employing purely repulsive WCA interactions between uncharged monomers. Lastly, the reaction MC moves in the GRM simulations efficiently insert/remove particles throughout the brush phase, facilitating free ion penetration and monomer ionization deep within the PEBs (Fig. S2–S5 in SI). As a result, the MC moves allow the simulations to bypass kinetic bottlenecks associated with chain ionization and ion diffusion inside the PEBs, which may prevent them from becoming trapped in long-lived metastable states like those observed in experiments. These three possibilities are not mutually exclusive, and we anticipate that using simulation to investigate the nature of swelling hysteresis may require developing explicit solvent PEB models with tunable hydrophobic interactions and simulation schemes that more closely mimic the dynamics of the experimental protocols.

Fig. 4 Brush height 〈h〉 during increasing (3 → 11) and decreasing (11 → 3) pH sweeps for a PEB with f = 1.0 at c_salt = 0.01 M. The sweeps were performed in an iterative fashion by using the final configuration from the simulation at pH^res_i to initialize the next simulation at pH^res_i+1. Uncertainties are smaller than the symbol size.

Examination of the brush height as a function of c^res_salt reveals that salt concentration also strongly affects the swelling behavior of the PEBs (Fig. 5). At pH^res = 11, the PEBs are fully deprotonated across the range of studied salt concentrations. As a result, they behave like neutral brushes and remain in relatively compact states as c^res_salt is varied. At pH^res = 3, the brushes are fully ionized and behave like strong PEBs, exhibiting two distinct regimes. In the osmotic regime at low c^res_salt, the PEBs are fully extended and the brush 〈h〉 exhibits little variation with c^res_salt. In the salted regime at high c^res_salt, by contrast, 〈h〉 begins to decrease monotonically with increasing salt concentration. This behavior reflects the collapse of the brush as electrostatic repulsions along the chains become strongly screened by the salt. The intersection of fits to the two limiting regimes allows for the crossover salt concentration to be estimated (Fig. 6). As expected, increases with f, reflecting the fact that higher salt concentrations are required to screen more highly charged PEBs. These trends are qualitatively similar to those observed in recent experiments on ionizable brushes (Fig. S6 in SI).⁴⁵

Fig. 5 Brush height 〈h〉 vs. c^res_salt at different pH^res for PEBs with ionizable monomer fractions f = (a) 0.2, (b) 0.6, and (c) 1.0. Uncertainties are smaller than the symbol size.

Fig. 6 (a) Brush height 〈h〉 vs. the reservoir salt concentration c^res_salt at pH^res = 3 for PEBs with different ionizable monomer fractions f. Dotted and dashed lines show fits to the osmotic and salted regimes, respectively. Uncertainties are smaller than the symbol size. (b) Osmotic-to-salted crossover concentration vs. ionizable monomer fraction f.

According to mean-field theory predictions for weak and strong PEBs,^10,11 the brush height in the salted regime should exhibit power-law scaling with the salt concentration 〈h〉 ∼ c^γ_salt with exponent γ = −1/3. We test this prediction by fitting 〈h〉 vs. c^res_salt to a power law in the salted regime at pH^res = 3, where the PEBs are fully charged. In each case, the observed scaling is weaker than predicted by theory (Fig. 6(a)). The theoretical scaling prediction of γ = −1/3 is derived based on the assumption that the concentration of salt is much lager than that of neutralizing counterions, such that the osmotic pressure is dominated by contributions from the former. Thus, the weaker than expected scaling behavior suggests that this limit has not been realized in the simulations. The maximum c^res_salt considered in the present study is ca. 1.0 M, beyond which the model fails to qualitatively predict the behavior of NaCl's activity in water.⁴⁸

Following ref. 26 and 27, we also examine the scaling of 〈h〉 with the total free ion concentration inside the brush c^brush_ion (Fig. 7). The species concentrations were calculated from the number of free ions in the volume L_x × L_y × z*, where z* is the upper bound on the brush phase along the z axis of the simulation cell, defined by the location of the inflection point in the monomer density profile. We observed that c^brush_ion > c^res_salt for all conditions examined (Fig. 7(a)). The extracted scaling exponent γ decreases with the ionizable monomer fraction f, approaching the predicted value γ = −1/3 as f → 1 (Fig. 7(b) and (c)). Thus, the high salt scaling regime is more easily reached by considering the total free ion concentration c^brush_ion rather than the salt concentration in the reservoir c^res_salt. Similar agreement with the theoretical scaling prediction of γ = −1/3 has been reported in experiments of weak PEBs under conditions where the monomers are fully charged²⁵ and in previous simulation studies of strong PEBs with f ≈ 1 in which the ion concentration inside the brush phase was used as the scaling variable.^26,27 Further analysis reveals, however, that even when using c^brush_ion as the scaling variable, the expected value of γ = −1/3 is only approximately recovered for pH^res ≤ 5 as f → 1 (Fig. 7(c)). Outside of these conditions, we find that |γ| < 1/3, indicating that the variation of brush height with salt concentration is weaker than predicted by mean-field theory. Weaker scaling relationships in the salted brush regime have also been observed experimentally and were attributed to excluded volume interactions due to monomer hydrophobicity and size.^24,45

Fig. 7 (a) Total free ion concentration inside the brush c^brush_ionvs. the reservoir salt concentration c^res_salt at pH^res = 3. (b) Scaling of the brush height 〈h〉 with c^brush_ion in the salted regime. Dotted lines indicate power-law fits to 〈h〉 ∼ c^γ_salt to estimate the scaling exponent γ. Uncertainties in (a) and (b) are smaller than the symbol size. (c) Scaling exponent γ vs. the ionizable monomer fraction f for different pH^res. The dashed line indicates the γ = −1/3 scaling exponent predicted by mean-field theory.

The brushes behave like strong PEBs and neutral brushes at pH^res = 3 and 11, respectively. However, at intermediate pH^res, the behavior is more complex. Specifically, in the osmotic regime at low c^res_salt, the brush height increases with salt concentration, reaching a maximum at before crossing over to the salted regime, where it again decreases with increasing c^res_salt. This non-monotonic behavior is characteristic of weak PEBs and arises from changes in the strength of the charge–charge repulsions along the chains resulting from competing ionization and electrostatic screening effects. We quantify the latter by computing the Debye length inside the brush via, where c^brush_i denotes the concentrations of charged species in the brush phase, which include both free ions and ionized monomers.

In the osmotic regime at low salt concentrations, increasing c^res_salt leads to an increase in brush ionization (Fig. 8(a)). For weakly acidic PEBs, it has been posited that added salt replaces protons within the brush and thus promotes ionization by shifting the local pH.⁷⁴ Analogous replacement of hydroxyl ions by salt may occur in weakly basic PEBs. Analysis of ion density profiles reveals enhanced anion uptake in the PEBs upon increasing c^res_salt (Fig. S4 in SI), but the effects of free chloride and hydroxide ions cannot be analyzed separately using the model because the two species have identical interactions with the polymer chains. Adding salt also increases electrostatic screening, as indicated by a decrease in the Debye length (Fig. 8(b)), which reduces the energetic penalty for ionization. Nonetheless, c^res_salt is not large enough in this regime to strongly screen electrostatic repulsions along the chain. Hence, the enhanced ionization effect is dominant and results in stronger electrostatic repulsions that stretch the chains as c^res_salt increases. The repulsions and chain stretching continue to increase with c^res_salt until they reach a maximum in the osmotic-to-salted crossover region. As c^res_salt increases further, the screening effect of the salt becomes dominant, diminishing the effects of charge–charge repulsions and leading to a concomitant decrease in brush height as c^res_salt increases.

Fig. 8 (a) Degree of ionization and (b) Debye length vs. salt concentration. Uncertainties are smaller than the symbol size.

The non-monotonic variation of brush height with c^res_salt has been observed in experiments.^23,45 It is also been reported by mean-field theory studies,^10,11 which predict that the brush height of weak PEBs in the osmotic regime should exhibit a power-law increase with the salt concentration 〈h〉 ∼ c^ν_salt with exponent ν = 1/3. Similar to the case for the salted regime (Fig. 7), we find weaker than expected scaling when c^res_salt is used as the scaling variable (Fig. 9(a)). The scaling exponent ν attains its largest values near pH^res = 6, where it increases from 0.02 to 0.06 as f increases from 0.2 to 1.0. These values are in reasonable agreement with the ν = 0.01–0.04 estimated from recent experiments.⁴⁵ The experiments also report the strongest scaling at intermediate pH 6–7, which is expected because the brushes behave like strong PEBs and neutral brushes at the pH extremes (11 and 3, respectively) and thus exhibit little variation in brush height with c^res_salt at low salt concentrations. By contrast, when the total free ion concentration inside the brush phase c^brush_ion is used as the scaling variable, we observe that the expected ν = 1/3 is approximately recovered at pH 6 as f → 1 (Fig. 9(b) and (c)). Hence, we find that the predicted scaling behavior in both the osmotic (ν = 1/3) and salted (γ = −1/3) regimes is recovered for weak PEBs with f = 1 when c^brush_ion is used in place of c^res_salt under appropriately chosen pH conditions.

Fig. 9 (a) Brush height 〈h〉 vs. the reservoir salt concentration c^res_salt at pH^res = 6 for PEBs with different ionizable monomer fractions f. (b) Scaling of the brush height 〈h〉 with c^brush_ion in the osmotic regime. Dotted lines in (a) and (b) indicate power-law fits to estimate the scaling exponent ν. Uncertainties in (a) and (b) are smaller than the symbol size. (c) Scaling exponent ν from power-law fits using c^brush_ion as the scaling variable (e.g., as in panel (b)) vs. the ionizable monomer fraction f for different pH^res. The dashed line indicates the ν = 1/3 scaling exponent predicted by mean-field theory.

4 Conclusions

We performed simulations using the grand reaction method (GRM) to investigate the responses of weak PEBs, consisting of polymer chains with random sequences of neutral and weakly basic monomers, to changes in solution conditions. The simulated titration curves decreased sigmoidally with increasing pH, qualitatively consistent with the behavior predicted by the Henderson–Hasselbalch (HH) equation. However, increasing quantitative deviations from the HH predictions were observed upon increasing the fraction of ionizable (weakly basic) monomers on each chain f, reflecting non-ideal behavior associated with charge–charge repulsions. Increasingly non-ideal behavior was also observed upon decreasing the concentration of monovalent salt due to a concomitant reduction in the screening of the electrostatic repulsions by free solution ions. This latter observation is consistent with experimental studies reporting a decrease in the apparent pK_a with decreasing salt concentration.^45,56,64

The brush height was also observed to decrease in a monotonic, sigmoidal fashion with increasing pH. The degree of chain stretching at low pH increased with f due to increased charge–charge repulsions, whereas the onset of swelling behavior shifted to higher pH as the salt concentration was increased due to the enhanced ionization resulting from electrostatic screening by free solution ions. The brush height response curves obtained upon incrementally increasing or decreasing pH were found to be indistinguishable. This reversible behavior is in sharp contrast with the pronounced hysteresis often observed in experiments, which has been attributed to the presence of long-lived metastable states arising from water-mediated charge-conformation coupling and transport limitations across the interface and within the brush. We posit that the absence of hysteresis in the simulations may be due to the PEB model's implicit representation of water and neglect of hydrophobic interactions, and the ability of GRM sampling protocol to efficiently equilibrate free ion and polymer charge distributions deep within the brush and thus circumvent kinetic bottlenecks associated with these metastable states.

Lastly, we found that the brush height exhibits a crossover from the osmotic brush regime to the salted brush regime, as predicted by mean-field theory and observed in experiments. As expected, the crossover solution salt concentration was found to increase with f due to the fact that higher concentrations of free ions are required to screen electrostatic repulsions within more highly charged PEBs. However, in both regimes, the scaling of the brush height with the solution salt concentration was observed to be weaker than predicted by mean-field theory. By contrast, when the measured salt concentration inside the brush was used as the scaling variable instead of the solution salt concentration, the theoretically predicted scaling behavior was approximately recovered in both regimes under specific conditions. The expected salted-brush scaling exponent was approached as f → 1 at pH ≲ 5, whereas the predicted osmotic-brush scaling exponent was recovered as f → 1 and pH → 6. Thus, agreement with mean-field theory was only observed using the brush salt concentration as the scaling variable for PEBs with f ≈ 1 at pH values below the intrinsic dissociation constant of the weakly basic monomer (pK^base_A = 8.4), where the chains are strongly ionized.

Our study highlights several outstanding challenges in understanding the behaviors of weak PEBs. First, experimental studies have used a variety of polymers. Differences in monomer chemistry and molecular weight have been posited to explain why some studies report scaling behavior consistent with mean-field theory, whereas others exhibit weaker than expected scaling behavior. Monomer chemistry effects are challenging to capture with the typical coarse-grained models used in PEB simulations studies. Second, different physical quantities are measured in experiment and simulation. For example, the brush swelling ratio (i.e., the extension of the chains relative to that in the PEB's dry state) is typically measured in experiments. In simulations, the absolute brush height is measured, and a dry reference state is not well defined in the context of the GRM. Similarly, agreement with mean-field theory has only been observed in simulations of strong and weak PEBs when the salt concentration within the brush phase is used as the scaling variable, but this quantity is not directly accessible in experiments. Lastly, the observed brush height hysteresis in many experimental studies indicates that these measurements are strongly influenced by presence of long-lived metastable states. Modified sampling protocols may be required to probe these metastable states in simulation to facilitate direct comparison with experiments.

Although the challenge in reconciling the observations from experiment, simulation, and theory currently hinder complete understanding of weak PEBs, they provide opportunities for future studies to make significant advances in this area. The effects of salt valency on PEB behavior is also an interesting avenue for future study. In the current work, we focused solely on symmetric monovalent salts, but experimental and computational studies have shown that salts with trivalent cations can induce collapse and structural inhomogeneities in strong PEBs.^75,76 A recent GRM study has also shown that divalent ions enhance the ionization of weak PEBs and lead to a two-stage pH swelling response.³⁶ Thus, GRM employed here provides an opportunity to investigate how multivalent and asymmetric salts influence charge regulation and pH effects in weak PEBs. Lastly, we have focused on examining “macroscopic” static equilibrium properties such as titration curves and the brush height response to facilitate comparison with theories for predicting these behaviors. Future studies of polymer and ion dynamics and how they couple to local heterogeneity within the PEBs (e.g., polymer chain and free ion distributions) are also of significant interest as they may provide insight into the origin of the kinetic bottlenecks associated with the pH response hysteresis observed in experiments.

Conflicts of interest

There are no conflicts to declare.

Data availability

Supplementary information (SI) is available including brush monomer, chain ionization, and free ion profiles and a comparison of the simulated brush swelling behavior with experimental measurements from ref. 45. See DOI: https://doi.org/10.1039/d5sm00999e.

The ESPResSo code used to perform the GRM simulations can be found at https://espressomd.org/.⁴⁶ The version of the code employed for this study is version 4.2.

Acknowledgements

This work was supported by the Welch Foundation (Grants E-1882 and the Welch Center for Advanced Bioactive Materials Crystallization, Award V-E-0001), the National Science Foundation (CBET-2113767 and DMR-2338550), and the ACS Petroleum Research Fund (New Directions Grant 67079-ND7). Computational resources were provided by the Hewlett Packard Enterprise Data Science Institute at the University of Houston and the Texas Advanced Computing Center at the University of Texas at Austin.

Notes and references

1. W. T. Huck, Mater. Today, 2008, 11, 24–32.
2. I. Tokarev, I. Tokareva and S. Minko, ACS Appl. Mater. Interfaces, 2011, 3, 143–146.
3. O. Wiarachai, T. Vilaivan, Y. Iwasaki and V. P. Hoven, Langmuir, 2016, 32, 1184–1194.
4. A. de los Santos Pereira, T. Riedel, E. Brynda and C. Rodriguez-Emmenegger, Sens. Actuators, B, 2014, 202, 1313–1321.
5. H. Vaisocherová-Lísalová, F. Surman, I. Víšová, M. Vala, T. Špringer, M. L. Ermini, H. Šípová, P. Šedivák, M. Houska and T. Riedel, et al. , Anal. Chem., 2016, 88, 10533–10539.
6. L. Zhang, H. P. Bei, Y. Piao, Y. Wang, M. Yang and X. Zhao, ChemPhysChem, 2018, 19, 1956–1964.
7. A. Kusumo, L. Bombalski, Q. Lin, K. Matyjaszewski, J. W. Schneider and R. D. Tilton, Langmuir, 2007, 23, 4448–4454.
8. G. Ferrand-Drake del Castillo, M. Koenig, M. Müller, K.-J. Eichhorn, M. Stamm, P. Uhlmann and A. Dahlin, Langmuir, 2019, 35, 3479–3489.
9. H. Gustafsson, A. Küchler, K. Holmberg and P. Walde, J. Mater. Chem. B, 2015, 3, 6174–6184.
10. P. Pincus, Macromolecules, 1991, 24, 2912–2919.
11. E. Zhulina, T. Birshtein and O. Borisov, Macromolecules, 1995, 28, 1491–1499.
12. Y. Tran, P. Auroy and L.-T. Lee, Macromolecules, 1999, 32, 8952–8964.
13. M. Balastre, F. Li, P. Schorr, J. Yang, J. W. Mays and M. V. Tirrell, Macromolecules, 2002, 35, 9480–9486.
14. H. Ahrens, S. Förster and C. A. Helm, Phys. Rev. Lett., 1998, 81, 4172–4175.
15. M. Biesalski and J. Rühe, Macromolecules, 2002, 35, 499–507.
16. M. Biesalski and J. Rühe, Macromolecules, 2004, 37, 2196–2202.
17. M. Biesalski, D. Johannsmann and J. Rühe, J. Chem. Phys., 2002, 117, 4988–4994.
18. T. Wu, P. Gong, I. Szleifer, P. Vlček, V. Šubr and J. Genzer, Macromolecules, 2007, 40, 8756–8764.
19. B. Lego, W. Skene and S. Giasson, Macromolecules, 2010, 43, 4384–4393.
20. H. Zhang and J. Rühe, Macromolecules, 2005, 38, 4855–4860.
21. J. D. Willott, B. A. Humphreys, G. B. Webber, E. J. Wanless and W. M. De Vos, Langmuir, 2019, 35, 2709–2718.
22. J. D. Willott, T. J. Murdoch, B. A. Humphreys, S. Edmondson, E. J. Wanless and G. B. Webber, Langmuir, 2015, 31, 3707–3717.
23. J. D. Willott, T. J. Murdoch, B. A. Humphreys, S. Edmondson, G. B. Webber and E. J. Wanless, Langmuir, 2014, 30, 1827–1836.
24. K. Ehtiati, S. Z. Moghaddam, A. E. Daugaard and E. Thormann, Macromolecules, 2021, 54, 3388–3394.
25. N. R. Hollingsworth, S. I. Wilkanowicz and R. G. Larson, Soft Matter, 2019, 15, 7838–7851.
26. C. Ibergay, P. Malfreyt and D. J. Tildesley, J. Phys. Chem. B, 2010, 114, 7274–7285.
27. N. A. Kumar and C. Seidel, Macromolecules, 2005, 38, 9341–9350.
28. C. E. Reed and W. F. Reed, J. Chem. Phys., 1992, 96, 1609–1620.
29. M. Ullner, B. Jönsson, B. Söderberg and C. Peterson, J. Chem. Phys., 1996, 104, 3048–3057.
30. M. Ullner and C. E. Woodward, Macromolecules, 2000, 33, 7144–7156.
31. A. Z. Panagiotopoulos, J. Phys.: Condens. Matter, 2009, 21, 424113.
32. F. Uhlík, P. Košovan, Z. Limpouchová, K. Procházka, O. V. Borisov and F. A. M. Leermakers, Macromolecules, 2014, 47, 4004–4016.
33. J. Landsgesell, C. Holm and J. Smiatek, J. Chem. Theory Comput., 2017, 13, 852–862.
34. S. Barr and A. Panagiotopoulos, J. Chem. Phys., 2012, 137, 144704.
35. K. Radhakrishnan, D. Beyer and C. Holm, Macromolecules, 2025, 58, 389–402.
36. D. Beyer and C. Holm, ACS Macro Lett., 2024, 13, 1185–1191.
37. D. Beyer, P. Košovan and C. Holm, Phys. Rev. Lett., 2023, 131, 168101.
38. J. Yuan and Y. Wang, J. Phys. Chem. B, 2021, 125, 10589–10596.
39. D. Beyer and C. Holm, J. Chem. Phys., 2023, 159, 014905.
40. D. Beyer, C. Holm and Z.-G. Wang, J. Phys. Chem. Lett., 2025, 16, 8245–8251.
41. X. Yuan, H. W. Hatch, J. C. Conrad, A. B. Marciel and J. C. Palmer, Soft Matter, 2023, 19, 4333–4344.
42. O. V. Rud, J. Landsgesell, C. Holm and P. Kosovan, Desalination, 2021, 506, 114995.
43. D. Beyer, P. Kosovan and C. Holm, Macromolecules, 2022, 55, 10751–10760.
44. J. Landsgesell, D. Beyer, P. Hebbeker, P. Košovan and C. Holm, Macromolecules, 2022, 55, 3176–3188.
45. F. Safi Samghabadi, S. Ramezani Bajgiran, M. Villegas Orellana, J. C. Conrad and A. B. Marciel, ACS Macro Lett., 2024, 13, 1570–1576.
46. F. Weik, R. Weeber, K. Szuttor, K. Breitsprecher, J. de Graaf, M. Kuron, J. Landsgesell, H. Menke, D. Sean and C. Holm, Eur. Phys. J.:Spec. Top., 2019, 227, 1789–1816.
47. R. Akter, X. Yuan, A. B. Marciel, J. C. Conrad and J. C. Palmer, J. Phys. Chem. B, 2025, 129(39), 10174–10183.
48. J. Landsgesell, P. Hebbeker, O. Rud, R. Lunkad, P. Kosovan and C. Holm, Macromolecules, 2020, 53, 3007–3020.
49. R. Stano, P. Kosovan, A. Tagliabue and C. Holm, Macromolecules, 2021, 54, 4769–4781.
50. K. Kremer and G. S. Grest, J. Chem. Phys., 1990, 92, 5057–5086.
51. R. B. Bird, R. C. Armstrong and O. Hassager, Dynamics of polymeric liquids, Wiley, 1987, vol. 1–2.
52. J. D. Weeks, D. Chandler and H. C. Andersen, J. Chem. Phys., 1971, 54, 5237–5247.
53. R. W. Hockney and J. W. Eastwood, Computer Simulation using Particles, 1988, pp. 267–304.
54. A. Arnold, J. de Joannis and C. Holm, J. Chem. Phys., 2002, 117, 2496–2502.
55. J. de Joannis, A. Arnold and C. Holm, J. Chem. Phys., 2002, 117, 2503–2512.
56. S. Ramezani Bajgiran, F. Safi Samghabadi, S. Li, J. C. Conrad and A. B. Marciel, Macromolecules, 2023, 56, 9218–9228.
57. B. Widom, J. Chem. Phys., 1963, 39, 2808–2812.
58. J. Jackson and L. Klein, Phys. Fluids, 1964, 7, 228–231.
59. W. Smith and B. Triska, J. Chem. Phys., 1994, 100, 3019–3027.
60. J. K. Johnson, A. Z. Panagiotopoulos and K. E. Gubbins, Mol. Phys., 1994, 81, 717–733.
61. V. S. Rathee, H. Sidky, B. J. Sikora and J. K. Whitmer, J. Am. Chem. Soc., 2018, 140, 15319–15328.
62. V. S. Rathee, B. J. Sikora, H. Sidky and J. K. Whitmer, Mater. Res. Express, 2018, 5, 014010.
63. V. S. Rathee, H. Sidky, B. J. Sikora and J. K. Whitmer, Polymers, 2019, 11, 183.
64. G. Ferrand-Drake del Castillo, R. L. Hailes and A. Dahlin, J. Phys. Chem. Lett., 2020, 11, 5212–5218.
65. D. Aulich, O. Hoy, I. Luzinov, M. Brücher, R. Hergenröder, E. Bittrich, K.-J. Eichhorn, P. Uhlmann, M. Stamm and N. Esser, et al. , Langmuir, 2010, 26, 12926–12932.
66. K. E. Secrist and A. J. Nolte, Macromolecules, 2011, 44, 2859–2865.
67. B. T. Cheesman, E. G. Smith, T. J. Murdoch, C. Guibert, G. B. Webber, S. Edmondson and E. J. Wanless, Phys. Chem. Chem. Phys., 2013, 15, 14502–14510.
68. J. Leyte, H. Arbouw-van Der Veen and L. Zuiderweg, J. Phys. Chem., 1972, 76, 2559–2561.
69. M. Ghasemi and R. G. Larson, Prog. Polym. Sci., 2021, 112, 101322.
70. V. Yadav, A. V. Harkin, M. L. Robertson and J. C. Conrad, Soft Matter, 2016, 12, 3589–3599.
71. N. Hollingsworth and R. G. Larson, Polymers, 2021, 13, 812.
72. S. C. Howard, V. Craig, P. A. FitzGerald and E. J. Wanless, Langmuir, 2010, 26, 14615–14623.
73. J. D. Willott, B. A. Humphreys, T. J. Murdoch, S. Edmondson, G. B. Webber and E. J. Wanless, Phys. Chem. Chem. Phys., 2015, 17, 3880–3890.
74. E. Zhulina and O. Borisov, Langmuir, 2011, 27, 10615–10633.
75. J. Yu, N. E. Jackson, X. Xu, B. K. Brettmann, M. Ruths, J. J. de Pablo and M. Tirrell, Sci. Adv., 2017, 3, eaao1497.
76. N. E. Jackson, B. K. Brettmann, V. Vishwanath, M. Tirrell and J. J. de Pablo, ACS Macro Lett., 2017, 6, 155–160.

---