2660
J. Phys. Chem. 1984,88, 2660-2669
composition x,,, = 0.1 123 are put into eq 15, and the temperature T satisfying eq 15 is calculated by an iterative procedure. This temperature is the temperature of phase separation of the HA/D,O mixture under consideration. The calculation is repeated for mixtures of different compositions (different values of xoHA and Y)using the same value of xaB. The slight change of xapwith Y is neglected. Figure 6 shows the result of these calculations. The calculated curve is shifted by 0.135 K to higher temperatures to overlap with the experimental data. The experimental data are the same as that shown in Figure 1. The same calculation is repeated for the system D A / H 2 0 with x,,, = 0.1 123. The result is shown in Figure 7. Here, the calculated curve is shifted to higher temperatures by 0.205 K. The experimental data are the same as in Figure 2. The agreement between the form of the calculated and experimentally determined curves in Figures 6 and 7 is satisfactory. The temperature shifts of the calculated curves are assumed to reflect the influence of impurities. The properties of the coexisting surface of the system isobutyric acid/H20, D 2 0 appear to be typical for three-component systems in which H / D isotope-exchange reactions take place. This can be concluded from measurements of the temperature of phase separation as function of composition in the system phenol ( O H ) / D 2 0 with an upper critical point” and the system 2-but-
oxyethanol ( O H ) / D 2 0 with a lower critical point.’* In both systems the mixture of critical composition has a lower D 2 0 content than the mixture of maximal and minimal phase-separation temperatures, respectively.2’ For phenol (OH)/D20, T, = 78.445 “C, XOHP,~= 0.095, T,,,,, = 78.495 “C, and x ~ H ~ = 0.085. For 2-butoxyethanol (OH)/D20, T, = 42.385 “C,xoHBp = 0.061, and Tp,min = 42.350 OC. Tables of experimental data shown in Figures 1 and 2 and data leading to eq 9 and 10 can be found in ref 2 1.
Acknowledgment. We thank G . Rottger, P. Harnisch, I. Schafke, and Dr. L. Belkoura for helping us during different stages of this study. The financial support of the Deutsche Forschungsgemeinschaft is also gratefully acknowledged. Registry No. Hydrogen, 1333-74-0; isobutyric acid, 79-3 1-2. (17) J. Timmermans and G. Poppe, C.R. Hebd. Seances Acad. Sci., 201, 524 (1935). (18) C. M. Ellis, J . Chem. Educ., 44, 405 (1967). (19) J. Zinn-Justin in “Phase Transitions”, M. Levy, J. LeGuillou, and J. Zinn-Justin, Eds., Plenum Press, New York and London, 1982, p 349. (20) W. J. Green, J . Chem. Eng. Data, 24, 92 (1979). (21) P. Gansen, Doctoral Thesis, University of Koln, 1984.
Ionic Atmosphere of Rodlike Polyelectrolytes. A Hypernetted Chain Study Russell Bacquet and Peter J. Rossky*t Department of Chemistry, University of Texas at Austin, Austin, Texas 78712 (Received: August 12, 1983; In Final Form: November 16, 1983)
Numerical solutions to the hypernetted chain (HNC) integral equation have been obtained for a model system representing an infinitely dilute rodlike polyelectrolyte in an aqueous 1-1 electrolyte solution. Bulk salt concentrations, cST, of lo-], and M, and reduced polyion charges, 6, ranging from 0.5 to 5.0, have been studied. The distribution functions are analyzed in terms of various structural parameters and compared with the counterion condensation (CC) formalism of Manning and with available Poisson-Boltzman (PB) data. Although HNC and CC are in general qualitative agreement, several quantitative aspects of CC theory, including the special nature off = 1.0 and the csT independence of the Manning radius, are not reproduced by HNC.
I. Introduction It is widely appreciated that the behavior of polyelectrolytes in solution is a sensitive function of the ionic environment of the polymer. Of special importance is the influence of solution composition on nuclei acid conformation’ and on the mode of interaction of nucleic acids with other solution components such as proteinsa2 A detailed understanding of the spatial distribution of small ions in the vicinity of a polyelectrolyte is fundamental to a microscopic interpretation of such phenomena. However, experimental studies of these distributions3 are, at present, only capable of elucidating the general character of the structure and frequently require substantial assumptions for their interpretation! Complementary theoretical descriptions are therefore of great interest. Theoretical treatments have been based, with few exception^,^^^ on either the counterion condensation (CC) formulation, developed primarily by M a n n i ~ ~ g , ~or- ’the ~ Poisson-Boltzmann (PB) e q ~ a t i o n . ’ ~ ,These ’~ approaches vary in the type and detail of information provided. C C theory is developed from a free energy minimization using a simple two-state model. The resulting description has the substantial advantage that results for a variety of phenomena depend on only a few physical parameters and can be obtained without the need for complex calculations. PB theory Alfred P. Sloan Foundation Fellow.
0022-3654/84/2088-2660$01.50/0
provides a more detailed description in that the small-ion distribution functions are obtained. However, the distributions obtained from PB theory account only for the mean electrostatic field due to the small-ion environment and inherently neglect other, short-ranged, forces among the small ions. Theoretical analysis15 indicates that, at least for highly charged polyelectrolytes in dilute (1) H. Drew, T. Takano, S. Tanaka, K. Itakura, and R. E. Dickerson, Nature (London), 286, 567 (1980). (2) A. Rich, N. C. Seeman, and J. M. Rosenberg in “Nucleic acid-Protein Recognition”, H. J. Vogel, Ed., Academic Press, New York, 1977. (3) See, for example: H. Magelhat, P. Turr, P. Tivant, M. Chemla, R. Menez, and M. Drifford, Biopolymers, 18, 187 (1979); J. L. Leroy and M. Gubron, ibid., 16, 2429 (1977); R. W. Wilson and V. A. Bloomfield, Eiochemistry, 18,2192 (1979); P. C. Karenzi, B. Meurer, P. Spegt, and G. Weill, Biophys. Chem., 9, 181 (1979). (4) C. F. Anderson, M. T. Record, Jr., and P. A. Hart, Biophys. Chem., 7, 301 (1978). (5) E. Clementi and G. Corongiu, Biopolymers, 21, 763 (1982). (6) D. Bratko and V. Vlachy, Chem. Phys. Lett., 90, 434 (1982); M. LeBret and B. H. Zimm, Biopolymers, 23, 274 (1984). (7) G. S. Manning, J . Chem. Phys., 51, 924 (1969). (8) G. S . Manning, J . Chem. Phys., 51, 934 (1969). (9) G. S . Manning, J . Chem. Phys., 51, 3249 (1969). (10) G. S. Manning, Biophys. Chem., 7, 95 (1977). (11) G. S. Manning, Biophys. Chem., 9, 65 (1978). (12) G. S. Manning, Q. Reu. Biophys., 11, 179 (1978). (13) A. Katchalsky, Pure Appl. Chem., 26, 327 (1971). (14) D. Stigter, J. Colloid Interface Sci., 53, 296 (1975). (15) M. Fixman, J . Chem. Phys., 70, 4995 (1979).
0 1984 American Chemical Society
,
~
~
~
Ionic Atmosphere of Rodlike Polyelectrolytes salt solution, the resulting errors may be small. Nevertheless, considering the high ion concentrations predicted by both C C and PB in the immediate vicinity of the polymer, it is reasonable to be concerned about the accuracy of such approximations in this particular, critically important, region. A review of the C C and PB theories as applied to D N A has appeared,16 and we do not go into further detail here. In the work reported here we apply the hypernetted chain (HNC) integral equation techniq~e,’~J* a method which has been very successful in the study of ionic liquids and solutions.19 In this approach the short-range forces and the correlation in position among small ions are included, although not exactly. The more detailed, and presumably more accurate, description of the polyion’s ionic atmosphere provides an improved basis for the study of more complex phenomena and is useful immediately as a comparison for less detailed but more convenient methods. It is the latter which we focus on here. Related calculations have been carried out for the distribution of ions in contact with charged planar surface^.^^^^ Although we use a somewhat different approximation, as discussed later, the success of these calculations as compared to the Monte Carlo simulation results of VaUleau and co-workersZ4is encouraging. In particular, deviations from PB results are found in that case, and the integral equation approach is capable of reproducing these differences. It should, however, be noted that this integral equation approach can lead to inaccurate results for the interaction between pairs of similarly charged polyelectrolyte m o l e c ~ l e s ? ~We ~ ~avoid ~ the polyion-polyion interactions here by considering only polyelectrolyte a t infinite dilution. W e adopt a polymer model corresponding to that studied by most earlier workers,7-1z~z7-31 namely, an infinitely long uniformly charged rod, and describe the solvent only through its dielectric constant. Such a model has been amply demonstrated to be a valuable caricature of both the polymer solutionl2 and the small-ion, supporting-electrolyte, solution.32 Of primary emphasis in the present work is a detailed comparison of the H N C description of structure with that implicit in the C C formalism of Manning, and with available Poisson-Boltzmann results both obtained for the same model. This comparison includes an analysis of the ionic distributions both in the immediate vicinity of the polyion as well as at large distance. We begin section I1 by describing the simple model system used to represent the rodlike polyion in electrolyte solution. We then discuss the H N C integral equation and its application in the present context. Our numerical results are reported in section I11 and comparisons with C C and PB data are made. A summary of the main conclusions and some remarks on future work are contained in section IV. We supply two appendices, one on the renormalized form of the integral equations and the other dealing (16) C. F. Anderson and M. T. Record, Jr., Annu. Rev. Phys. Chem., 33, 191 (1982). (17) R. 0.Watts in ”Specialist Periodical Reports, Statistical Mechanics”, Vol. 1, The Chemical Society, London, 1973. (18) J. M. J. Van Leeuwen, J. Groeneveld, and J. DeBoer, Physica, 25, 792 (1 959). (19) D. N . Card and J. P. Valleau, J . Chem. Phys., 52, 6232 (1970); J. P. Valleau, L. K. Cohen, D. N. Card, ibid., 72, 5942 (1980). (20) D. Henderson, L. Blum, and W. R. Smith, Chem. Phys. Lett., 63, 381 (1979); D. Henderson and L. Blum, J. Electroanal. Chem., 111,217 (1980). (21) S.L. Carnie, D. Y.C. Chan, D. J. Mitchell, B. W. Ninham, J. Chem. Phys., 74, 1472 (1981). (22) M. Lozado-Cassou, R. Saavedra-Barrera, and D. Henderson, J . Chem. Phys., 77, 5150 (1982). (23) M. Lozado-Cassou and D. Henderson, J. Phys. Chem., 87, 2821 (1983). (24) G. M. Torrie and J. P. Valleau, J. Phys. Chem., 86,3251 (1982), and
references therein. (25) G. N . Patey, J . Chem. Phys., 72, 5763 (1980). (26) M. Teubner, J. Chem. Phys., 75, 1907 (1981). (27) H. P. Gregor and J. M. Gregor, J . Chem. Phys., 66, 1934 (1977). (28) J. A. Schellman and D. Stigter, Biopolymers, 16, 1415 (1977). (29) L. Kotin and M. Nagasawa, J. Chem. Phys., 36, 873 (1962). (30) D. Stigter, J . Phys. Chem., 82, 1603 (1978). (31) M. Gubon and G. Weisbuch, Biopolymers, 19, 353 (1980). (32) See, for example, P. S.Ramanathan and H. L. Friedman, J. Chem. Phys., 54, 1086 (1971).
The Journal of Physical Chemistry, Vol. 88, No. 12, 1984 2661 with the electroneutrality moment condition for rod-ion correlation functions. 11. Methods A . Model. The polyelectrolyte model chosen for the present study consists of an isolated and infinitely long cylinder bearing a uniform axial charge density given in reduced units by
where b is the contour length incorporating one electronic unit of charge, Tis the temperature, kB is Boltzmann’s constant, e is the dielectric constant of the pure solvent, and e is the magnitude of the electronic charge. Throughout this study we have used t = 78.358 (characterizing water) and T = 298.16 K so that f = 7.1507/b (2) where b is in angstroms. f = 4.2 is the value representative of double-stranded DNA. With an isolated infinite cylinder, end effects and polyionpolyion interactions do not enter. A cylindrical model is appropriate for polyions consisting of locally rodlike segments which are straight on the scale of the Debye ionic screening l e r ~ g t h , ~ ” ~ ~ which is the case for DNA.lZ The continuous axial charge distribution is a first approximation expected to be adequatez9 except in the case of low polyion charge density, where the equivalent discrete charge distribution is necessarily sparse. This case is not, however, of interest here. In our solution model the solvent is treated as a dielectric continuum. There are two species of small spherical ions which differ here only in the sign of their charge and which have “soft” repulsive cores. The time-independent equilibrium solution structure is expressed in terms of a s~lvent-averaged,~~ pairwise-additive,potential. Between the small ions we use an effective pair potential of the form32
where the i and j subscripts denote either of the small ion species, Al is a chosen constant, Zi is the signed charge in units of e for species i, ri* is an effective diameter for species i, and rij is the distance between ionic centers. For the effective pair potential between a small ion and an infinite charged cylinder bearing a uniform negative charge we use
where the p subscript denotes the polyelectrolyte, and r is the radial distance in a cylindrical coordinate system centerefon the rod. Throughout this study we have chosen A , = 7.4245 X erg A Az = 1.0383 X zi = f l
erg A
A uniform dielectric constant, extending through the rod, is implicit in these choices for the interactions in the model solution. The short-range r4 potential given in eq 3 and 4 has been used p r e v i o u ~ l yin~ ~small-ion * ~ ~ electrolyte solution studies as a more realistic alternative to a hard-sphere potential. We have adopted a corresponding form for the rod/ion core repulsion, and no special significance should be attached to the particular choice of power law employed here. In previous polyion work,9J0J2J7-31an impenetrable, hardsphere-like, potential has typically been employed for the small (33) A. Katchalsky, Z. Alexandrowicz, and 0. Kedem in “Chemical Physics of Ionic Solutions”, B. E. Conway and R. G. Barradas, Ed., Wiley, New York, 1966. (34) L. Onsager, Ann. N.Y. Acad. Sci., 51, 627 (1949). (35) W. G. McMillan and J. E. Mayer, J . Chem. Phys., 13, 276 (1945). (36) P. J. Rossky, J. B. Dudowicz, B. L. Tembe, and H. L. Friedman, J . Chem. Phys., 73, 3372 (1980).
2662
The Journal of Physical Chemistry, Vol. 88, No. 12, I984
ion-polymer interaction. However, it has been shown36for corresponding simple electrolyte solutions that there is a correspondence in both thermodynamics and solution structure between hard- and soft-sphere models if they are matched appropriately. The matching requirement is that the minimum in the pair potential for oppositely charged species occur at the same distance for both systems. We have used this same criterion for matching hard and soft models of the sphere/rod system, although no studies of the correspondence have been done for this case. The soft model system chosen in the present study was designed to correspond to a hard sphere/hard rod system in which the hard sphere had a radius of 1.38 8, (intermediate between Na+ and Cl-) and the center-to-center rod-ion contact separation was 10.0 8, (roughly representative of doubly stranded DNA). These considerations led to the following parameters for the “soft” model system
rp* = 8.909 8,
dependence by virtue of the infinite extent and uniformity of the rod. Denoting by kr the radial cylindrical component in k-space, we have immediately
Lpi(kr) = Epi(kr)
+ Fpj?pj(kr) Rij(kr)
(8)
which is the form of the O Z equation used in the present calculations. The Fourier transforms are carried out by using a previously described technique designed specifically for very long-ranged correlation functions.38 The usual measure of accuracy,38the electroneutrality defect eq B13 in Appendix B), was less ifi all cases, indicating that quite accurate numerical than solutions to the integral equations were obtained. The H N C integral equation consists of eq 5 together with the usual approximate closure relation39 (in coordinate space)
+
ri* = 0.934 8, One feature of the chosen pair potentials (eq 3 and 4) which should be noted here is that the strength of the short-range repulsion depends on the charges of the species. The (desired) result of this dependence is that, with all other parameters fixed, the minimum in the total pair potential for oppositely charged species occurs a t the same distance regardless of the charge magnitudes involved. The rod becomes “softer” as its charge is decreased and the Coulomb attraction weakens. This change in the core potential is not chosen for realism but rather to allow a closer comparison with earlier hard-sphere-model studies, while retaining some of the more realistic modeling associated with a soft core. The position of the first peak in the counterion distribution would shift somewhat to decreasing r with increasing Coulomb charge were a fixed core potential employed. B. HNC Integral Equation. For an isotropic multicomponent system in which one species is infinitely dilute (in this case, the polyelectrolyte) the Ornstein-Zernike ( O Z ) equation3’ quite generally consists of the set of coupled integral equations
where the subscript x can represent any of the species, the subscripts i and j denote any of the species except the one at infinite dilution, p j is the number density of species j, and rx represents, in general, the set of all coordinatas needed to specify the position of a “molecule” (atom, molecule, rigid body, etc.) of species x. hxiand cXiare the well-known total pair correlation function and direct correlation function, respectively. If species x is one of the small ion species, then eq 5 is the equation for a pure system (no polyelectrolyte, or other infinitely dilute species, present). Thps, in the polyelectrolyte studies, the simple salt results are first obtained and then incorporated through the term hj, on the right-hand side of eq 5. In the present case, species x is the infinitely dilute rodlike polyelectrolyte. The small-ion positions are specified by a three vector in, for example, Cartesian coordinates, and the rod position is specified in addition by orientation parameters. To obtain the working equations for the calculations here, we can assume that the rod lies, for example, along the z axis of our coordinate system, passing through the origin. Thus, we have, in this coordinate system (rij = ri - rJ)
Using the convolution theorem, one has, equivalently 1
Bacquet and-Rossky
(7)
where we have indicated the Fourier transform of a function by a caret. Adopting cylindrical coordinates, we have no angular dependence in the correlation functions by symmetry, and no z (37) L.S.Ornstein and F. Zernike, Proc. Acad. Sci. Amsterdam, 17,793 (1914).
cXiHNC = hxiHNC - In (hxiHNC 1 ) - fluxi (9) where 0 is l / k , T . It is not necessary to use the approximate closure in eq 9 for all components x in the system. In general, one may employ any source for the bulk correlations h,i in eq 6 . In particular, Henderson and c o - w ~ r k e r shave ~ ~ , emphasized ~~ the alternative approximation, introduced by Carnie et al.,2’ which uses the mean spherical approximation (MSA) for the bulk correlation functions h,, but the H N C closure to relate cpIto hpl. Reasonable agreement has been demonstrated between the hybrid approximation and simulation result^^^^^^ for the distribution of small ions near a charged plane, while results for the same system using H N C bulk correlation functions have been considered only briefly.20 However, it is clear that the H N C description of the bulk correlations is quantitatively superior to the MSA,19 so that it appears necessary to attribute the success of the hybrid theory to a fortuitous cancellation of errors.22 Whether such cancellation would pertain in the present case is unclear. Hence, we have chosen to use the best available results for the bulk correlation functions, namely, those from the H N C approximation, and correspondingly use the closure expressed in eq 9 throughout. Finally, it is important for computational purposes, as well as for analysis, to recast the integral equations in a form in which the long-ranged bare Coulomb potential is absent. In its place the Debye-Huckel (DH) screened potential appears. This can be accomplished in a straightforward manner, aqalogous to the well-known corresponding rearrangement for simple atomic salt solutions. These “renormalizations” are outlined ih Appendix A. 111. Results Calculations have been carried out for stoichiometric salt concentrations (csT) of 0,001, 0.01, and 0.1 M. For each concentration, eight values of rod charge in the range from 5 = 0.5 to f = 5.0 were examined. In the following, we first consider the counterion distribution near the rod, with particular emphasis on the amplitude of the local concentration. We focus on the significance of both the rod charge and the salt concentration in determining this feature. We then consider the long-range behavior of the counterion distribution, again as a function of cSTand f , and compare the HNC results with Debye-Huckel-like approximations which employ an effective rod charge.12 Next we examine the co-ion density profiles, with special attention to the f-dependence at short range. Finally we analyze the H N C results in terms of structural parameters presented by other authors and compare with C C and PB data. A . Counterion Distributions. 1 . Short-Range Behavior. The representative counterion correlation functions hp+(r),for csT = 0.1 M and various are presented in Figure 1. hp+(r)changes smoothly as f is increased through unity, which is the “critical” value marking the onset of condensation in Manning’s theory.12 We note that the small variations observed near 8 A are due to the implicit change in strength of the core repulsion (eq 4) as (38) P. J. Rossky and H. L. Friedman, J. Chem. Phys., 72,5694(1980). (39) H.L.Friedman and W. D. T. Dale in “Statistical Mechanics”, Part A, B. J. Berne, Ed., Plenum Press, New York, 1977.
,
The Journal of Physical Chemistry, Vol. 88, No. 12, 1984 2663
Ionic Atmosphere of Rodlike Polyelectrolytes
2. 8 0
2 . 10
1 . 40
0.70
0.00 6 . 00
14. 00
22.00
00
30. 0 0
rc61
~ C R I
Figure 1. Rod-counterion pair correlation functions hp+(r)vs. radial rod-ion center-to-center distance, r, for cST = lo-’ M. The curves correspond, in order of increasing amplitude, to rod charges o f f = 0.5,0.9, 1.0, 1.1, 2.0, 3.0, 4.2, 5.0.
0.28
Figure 3. Local concentrations of counterions ~ + ~ ~ as( in r )Figure ; 2 but for 6 = 2.0, 3.0, 4.2, 5.0. TABLE I: Maximum in c+m(r) E 10-3 M 0.5 0.009 0.9 0.038 1.o 0.05 1 1.1 0.067 2.0 0.357 3.0 0.984 4.2 2.092 3.013 5.0
(M)
: 4 I
0.21
0. 1 4
c b-I 0.07
0. 00 6. 0 0
14.00
22.00
30. 0 0
rc61 Figure 2. Local concentrations of counterions ~ + ~ ~for (cSTr =) lo-’ (-), 10” (--) and lo-’ (---) M. For each value of csT, the four curves correspond to the values 6 = 0.5,0.9, 1.0, and 1.1 (concentration in molar units).
discussed in section IIA. The results for cST = 0.01 and 0.001 M have been found to be very similar in form to those in Figure 1 except for a scale factor equal to the salt density so that p+hp+ is roughly constant. This scaling behavior is an alternate reflection of the near cST independence of the local concentration of counterions near the rod. More precisely, we define the local concentration via
where we recall that p+ is the bulk counterion concentration, and p+hp+is the excess counterion concentration. c+Lw for all three bulk salt concentrations and various f are shown in Figures 2 and 3. In Figure 2, we have grouped those results for f 5 1, while in Figure 3 results for f 1 2 are shown. The asymptotic value of each curve is, of course, trivially predictable. However, in the region near the rod, the behavior changes qualitatively as the rod charge density increases. Although a strong CST dependence is observed a t the lower f values, it
lo-* M
lo-’ M
0.034 0.082 0.100 0.120 0.425 1.018 2.034 2.878
0.160 0.240 0.265 0.291 0.625 1.190 2.112 2.858
diminishes with increasing f . At the highest rod charges (f 2 4), the local counterion concentration in the vicinity of the rod changes by only a few percent for bulk concentration changes of two decades. In general, however, the value at the maximum in c+Lw increases with cST for fixed [ except at the highest rod charge examined, where in fact this trend is reversed although the differences are small. The maximum values of ~ + ~ obtained ~ ( r )are summarized in Table I. It is clear from these values that although the local concentrations are rather independent of cST at high E, the value f = 1 plays no special role in the transition to this behavior. The values increase in nearly parallel fashion, with the differences remaining nearly the same up to [ 4. We note here also that the values obtained in the region very hear r = 10 A (rod-ion “contact”) are significantly smaller than are those obtained with a corresponding hard-sphere-like interaction poter~tial.~’This fact is totally and results from the artificiality of the infinite discontinuity in this latter form. 2. Long-Range Behavior. The counterion distributions decay quickly to their bulk values (Le., &OC(r) = cST, hp+= 0 ) , a result expected from the Debye-Huckel theory. Of most interest, however, is the question of the precise form of this decay. In the Bjerrum-like40 theory of Manning, the counterions far from the polyion are presumed to obey Debye-Hiickel-like (DH) behavior, but with an effective polyion charge density given by teff= 1.‘J2 Representative results pertaining to this question are shown in Figures 4-6, where we plot in logarithmic form counterion correlation functions from the H N C theory as well as analytical Debye-Huckel-like approximations. We note at the outset that for the H N C results, the amplitude of the counterion correlation function hp+(r)is not, in general, monotonically increasing with increasing [ for large r (r R 30 A), an apparent result of the detailed treatment of ionic screening.
-
(40) N.
Bjerrum, K.Dun. Vidensk. Selsk., 7, No. 9 (1926).
2664 The Journal of Physical Chemistry, Vol. 88, No. 12, 1984 - 1 . 60
Bacquet and Rossky
1
k -2.
GO -2. 60
-3. G O
-3. 6 0
In h
In h
-4. G O
-4. 6 0 -5. G O
80.00
1 1 0 . 00
140.00
170.00
rCfl1 Figure 4. Comparison with Debye-Hiickel approximations for In hp+(r) for CST M; 5 = 0.5 (X), 0.9 (a), 1.0 (I),1.1 (*), 2.0 ( O ) , 3.0 (a), 4.2 (+), 5.0 ( 0 ) ;In hp+DH(r) for teff= 1.0 (---) and In hp+DHL(r) for
teff= 1.o (--).
t -0. 500
-1.000
In h -1. 50
-2. 00
I
I
-5. 60 80. 00
110. 00
1 0. 00
140. 00
rtil Figure 6. Variation of tefffor the Debye-Huckel approximation, In hp+DH(r),for cST = lom2M. The dashed lines correspond, in order of with teff= 0.5, 0.9, 1.0, 1.1, 2.0; solid increasing amplitude, to In hPDH(r) lines are In hp+(r)with 5 = 0.5 (X), 1.0 (E), and 4.2 (+). in Appendix A). In eq 11, a is the center-to-center contact distance of a hard rod and a hard ion and has been set to 10 A, the position of the counterion maximum in our corresponding soft model. As is clear from Figures 4 and 5, eq 11 and 12 differ very little for the present purposes. The nearly parallel curves in the figures show that the rate of decay is essentially the same, and DH-like, for all cases at a given salt concentration. It is the amplitude which is more sensitive to f . At lo-* M,the D H results, with teff= 1, lie between the M (where the H N C f = 1 and f = 2 results, while for amplitude of hp+is much larger at the same value of r), the D H results appear to reproduce the higher f curves well. It is important to examine the sensitivity of the D H result to the choice of teff,and this is shown in Figure 6 for the case cST = M. teff= 1.2 provides a good fit with the high-f H N C results, and a “best” value not far from unity, for both counterions and co-ions, appears appropriate for the other concentrations as well. Note that in Figure 6, the H N C result for f = 0.5 is matched well by the D H result based on teff= f . This is in general agreement with C C theory; Le., the lack of “condensation” and no special reduction of the rod’s effective charge for 4 < 1.0.’J2 It is, however, interesting to note that the 5 = 1 curve is not particularly well represented by the D H result at this lower concentration; the D H result manifests significantly less screening than is found with H N C . From the behavior shown in Figures 4-6 as a function of $, and teffit is clear that quantitative reproduction of the long-range correlations requires a variable fa. Nevertheless, the concept that, for large f , a D H theory with teff 1 is a reasonable approximation is certainly supported, although, as for the short-range structure, f = 1 does not play a special role except as a qualitative division between “large” and “small” f . B. Co-ion Distributions. The co-ion results are displayed in Figures 7-9. For all three concentrations, increases in the rod charge initially correspond with decreases in the co-ion density for all r values. After f reaches some intermediate, csrdependent, value, this trend is reversed for a large range of r values; Le., co-ion density increases with increasing f . The magnitude of the increases is largest a t the lowest concentrations. There is a particular enhancement of h, in a small region centered at 15 A. At f = 5.0 we find a distinct peak in this region for cST = M and a slight bulge even for cST= lo-’ M. We note that these features actually involve only a very small fraction of the total ion density as is evident from the small amplitudes of h,(r) (cf. Figure 1).
-
-
(41)
T.L.Hill, Arch. Biochem. Biophys., 57,229 (1955).
The Journal of Physical Chemistry, Vol. 88, No. 12, 1984 2665
Ionic Atmosphere of Rodlike Polyelectrolytes
-0.4 4 0
0.000
-0.250
-0.580
-0.5CO
~
-0. 720
hF1
h
-0.750
w1
-0. 860 -1.000 00
16. 0 0
26. O@
3-. 00
Figure 7. Polyion-co-ion radial correlation function h,(r) for cST= lo-’ M ( 4 values and symbols as given in Figure 4). -0.250
Figure 9. Polyion-co-ion radial correlation function h,(r) for cST= lo-.’ M. -2. 00
-0.500
hF I
-3. 00
-0. 7 5 0
-4.00 6.00
16. 00
3”. 0 0
2 6 . 00
In
h
-5. 00
rCR1 Figure 10. Comparison with Debye-Hiickel approximations for In (-h,(r)) for cST = lo-* M. Curves are indicated as in Figure 4.
shell containing the condensed ions is independent of salt concentration and polyion radius and, for univalent counterions, is given by10-12(angstrom units) Vp = 68.23([ - l)b3 (13) The Manning radius, RM‘, is the radial thickness of V, plus the radius of the rod. Using eq 13, we have RM’ = [21.72([ - l)b2 a2]1/2 (14)
+
We use RM (no prime) to mean the actual distance enclosing a charge of OM for a given theory (e.g., HNC or Poisson-Boltzmann), reserving the symbol RM’for the radius predicted by eq 14 of CC theory. We define the integrated charge from species i, qi, by qi(r) = Zipixbx2TJr[h,i(r’)
+ lIr’dr’d4
dz
(15)
and the total integrated charge by ~~
(42)R. Bacquet and P. J. Rossky, J . Chem. Phys., 79, 1419 (1983).
qTOT(r) = Cqi(r) I
(16)
2666 The Journal of Physical Chemistry, Vol. 88, No. 12, 1984 1.000
j
p
-
Bacquet and Rossky 1.000
1
I
1 0.750
-
0.7 5 0
0.500
i
0. 500
qb-I
1
q (0
0.250
-
0. 2 5 0
0.000 6 . 00
1 4 . 00
22. 03
0.0 0
3 0 . 00
100. 0 0
L1
200. 00
300. 0 0
rCFl1
r C i 1 Figure 11. Total integrated charge qmT(r)(eq 16) for cST = lo-' M. The curves in order of increasing amplitude (beyond 10 A) correspond t o t = 0.5, 0.9, 1.0, 1.1, 2.0, 3.0, 4.2, 5.0; the point qToT(r=RM)= 1 rl is marked (*) for each 5 > 1.0.
Figure 13. Total integrated charge qToT(r)for csr = lo-' M. Notation as in Figure 11. 1. 0 0 0
4(rl
! . 000 t
0. 5 0 0 C.750
0.000
0.500
1 q Irl
c r -
0. 500
1 -
3.000
0.00
4C. C0
EC.
CO
123. 30
r C i 1 Figure 12. Total integrated charge qToT(r)for cST = lo-* M. Notation as in Figure 11.
qTOT(r) is the total net ionic charge inside a cylindrical volume of radius r and height b centered on the rod. Since b is the distance along the rod which incorporates one unit of electronic charge, electroneutrality requires that
Figures 11-13 present H N C results for qT&) for a number of concentrations and rod charges. Each constant-5 curve is marked with an asterisk at that point where r = RM,i.e., qT@) = 6 ~ . The speed with which qTOT(r) approaches its limiting value of unity is clearly very much dependent on bulk salt concentration. There is also such a dependence in the relation between R M and 5, as can be seen in the figures. In lo-' M salt, R M always increases with increasing 5 for the range of rod charges studied. For the lower concentrations, R M increases to a maximum value and then declines slightly as 4 is increased. The sign of the change in R M with 4 is determined by the imbalance in two offsetting trends. An increase in rod charge leads to an increase in counterion density near the rod and also to an increase in the 6~ needed to reduce the rod's effective 5 to 1. It is clear that for any given set of parameters this balance may have either sign. It is of interest to note that RM in C C theory exhibits such a nonmonotonic behavior as well (see Table 11). The H N C results may be readily
0.000 0.00
28.00
r
28.00
0.00
~
i
56. 00
r~ ~ k
Figure 14. Total integrated charge, qToT(r),for cST= lo-' ( 0 ) , (0), and (A)M and = 1.1 (A), 2.0 (B), 3.0 (C), 4.2 (D); solid line is 1and dashed line is RM' from eq 14.
compared graphically with condensation-theory predictions of R M , as we do in Figure 14. The horizontal solid line is 6~ for the given f'. The constant-csT H N C curves, qTOT(r),intersect the 6~ line at R M and can be compared to RM'which is indicated by the vertical dashed line. RM' is calculated from eq 14 with a = 10 A. Contrary to condensation theory7-12 and in qualitative agreement with PB data,31the H N C results show that R M is cST dependent. The condensation-theory RM' is smaller than the H N C RM except at the lower rod charges. The R M and RM' values from Figure 14 are collected in Table 11, along with R M based on counterions only, q+(r),RM' based on an alternative definition of a,and R M from a numerical solution by Gueron and Weisbuch31 of the PB equation (for an impenetrable rod, a = 10 A). In part B of Table 11, we include the result of application of the CC-theory RM' values to the H N C results. This gives a set of enclosed charge fractions which we denote as BH and a set of corresponding effective rod charges teff. Considering first the definition of rod radius, we note that the RM' of condensation theory depends on the radius, a, of an impenetrable rod (although V does not). A comparison of condensation theory with the sod-rod HNC results requires a careful choice of a. In Figure 14, we obtained RM' by setting a = 10 A,
The Journal of Physical Chemistry, Vol. 88, No. 12, 1984 2667
Ionic Atmosphere of Rodlike Polyelectrolytes TABLE II: Structural Parameters
A. Manning Radius RM.A
t 1.1 1.1 1.1 2.0 2.0 2.0 3.0 3.0 3.0 4.2 4.2 4.2 5 .O 5.0 5.0
CST,
M
10-3 10-2
lo-’ 10-3 10-2 10-1
10-3 10-2
lo-‘ 10-3 10-2 lo-’ 10-3 10-2
lo-’
condensation theory CLOC (a) I 0 . k 13.85 11.8 13.85 11.9 12.0 13.85 18.1 19.43 18.2 19.43 18.3 19.43 17.3 18.62 17.4 18.62 17.5 18.62 17.36 16.1 17.36 16.1 17.36 16.2 16.66 15.3 16.66 15.4 16.66 15.5
PoissonBoltzmann“
a = 10 8,
43 16.5
60 21.5
HNC based on qTOT based on q+ 14.10 14.10 11.52 11.51 9.84 9.77 36.05 35.97 23.34 23.08 15.78 15.29 38.99 38.86 26.88 26.46 18.39 17.64 35.05 34.94 27.35 26.86 19.78 18.86 3 1.04 30.93 26.73 26.21 20.24 19.24
B. Charge Fractions and Effective Rod Charges OHb
1.1 2.0 3.0 4.2 5 .O
lo-’ M 0.35 0.64 0.67 0.70 0.7 1
OFrom ref 31. b 8 is~ defined by 8~ = 8 ~ then , teff= 1.0.
lo-’ M 0.16 0.43 0.54 0.63 0.67 qTOT(RM’) = 6H
[eft
10-3 M 0.09 0.35 0.50 0.63 0.68
OM
0.091 0.500 0.667 0.162 0.800
lo-’ M 0.71 0.72 0.99 1.26 1.45
IO-’ M 0.92 1.14 1.38 1.56 1.65
10-3 M 1 .oo 1.30 1.50 1.56 1.60
where RM’ is obtained from eq 14 with a = 10 A. ‘Obtained from teff= (1 - 8 ~ ) tif;
the position of the minimum in the rod-counterion potential. It might be argued that a hard rod of radius less than 10 8,would correspond more accurately with our soft rod. To test this, Table I1 lists the values of RM’ obtained by choosing a to be the distance at which cLoc(r) = O.lcsT. We find that this rather extreme definition for a reduces RM’ by only 1-2 8, and usually increases the disagreement between H N C and CC. Another potential source of modification arises from the fact that, in Manning’s theory, 6 M refers only to counterion density. Therefore, q+(r) may by more appropriate than qTOT(r)for obtaining RM. However, as shown in Table 11, there is very little difference between the two functions for r 2: RM since q-(r) is still very small a t this distance. Comparing the available Poisson-Boltzmann results3*with both C C and HNC, we see that in general PB and H N C appear reasonably close to one another. They, however, agree most closely where all three approaches are in general agreement. At the lowest concentration (cST= M) PB and H N C agree that RM >> RM’ but the PB values are substantially larger than those from HNC. At this stage, it is not possible to state which is more accurate, but the fact that H N C is a more detailed theory, and includes small ion correlations, would suggest more confidence in this result. Table 1I.B provides a convenient basis for a qualitative comparison of H N C and CC. As is clear from these values, a general agreement of charge fractions between the two approaches is apparent, particularly for higher f. Hence, despite its quantitative shortcoming, there is little question that the C C theory does capture the qualitative structural features.
IV. Summary This paper has reported an application of the hypernetted chain integral equation to a simple model system representing an isolated rodlike polyion in a 1-1 supporting electrolyte solution. Rod-ion radial distribution functions have been obtained and analyzed for several bulk salt concentrations and a wide range of polyion charges. The results for the counterion distribution near the rod exhibit the strong direct dependence on f and insensitivity, except at low
f, to bulk salt concentration which has long been predicted theoretically and observed experimentally. The counterion distributions predicted by H N C change in a continuous fashion as f increases through unity and are always smooth functions of distance. Correspondingly, we find no evidence that a special character should be assigned to the value f = 1 .O. N o obvious basis for the physical division of counterion density into “condensed” and “free” portions for f > 1.0 is apparent, although we have not pursued this avenue in much detail here. A two-state theory consistent with the PB equation has been derivedlS and a similar possibility exists for the H N C formulation. In accord with C C theory,’J2 the long-range portions of h,+(r) and hp-(r)are well approximated, particularly at low salt concentration, by a D H approximation in which the polyion is assigned an effective reduced charge of either 4‘ (the physical charge) or unity, whichever is smaller. This approximation is not quantitative, however. Since RM encloses the quantity of charge necessary to achieve teff= 1.O,and the usefulness of this condition is well established, we have compared the values obtained by the various methods. Except at low l , the R M ’ of the CC theory,’0-12which is independent of salt concentration, is low compared to the H N C RM values, and the latter manifest significant cSTdependence. The PB results3’ for RM agree well with H N C at cST= lo-’ M but are substantially larger than H N C at M. Although a definitive test of the accuracy of the various results must await further study, this difference between H N C and PB can be reasonably attributed to the inclusion of small-ion correlations in the former approach, and the correspondingly higher counterion density (ionic atmosphere contraction) which would be expected to result. In summary, the HNC results suggest that the structural picture invoked in C C theory has many qualitatively accurate features, particularly at the high charge densities characteristic of nucleic acids. However, a number of quantitative aspects of the theory including the special nature of f = 1, and the constancy of RM with salt concentration, are not borne out here. This suggests that the successes of the theory may in part be a result of the focus on the high-f experimental examples which are of most interest.
2668 The Journal of Physical Chemistry, Vol. 88, No. 12, 1984
Bacquet and Rossky
Currently, computer simulations for selected cases are under way in this laboratory in order to test the quantitative accuracy of the present approach. At the same time, extensions to the important case of salt solutions composed of a mixture of ions of various valence are being carried out and will be reported in a future publication.
Q can be easily computed from eq A9 by using the Fourier transform, Fourier-Bessel transform, and the convolution theorem. q’can be obtained from eq A7 and has been shown9 to have the form
Acknowledgment. Support of the research reported here by a grant from the National Institutes of Health is gratefully acknowledged. P.J.R. is the recipient of a Research Career Development Award (PHS CA00899) awarded by the National Cancer Institute, DHHS.
where KOis the zeroth-order modified Bessel function of the second kind. When we compare eq A1 1 with eq 12, it is clear that the Debye-Huckel approximation for h ’is
Appendix A Renormalized Form of Integral Equations. For the case of a pure atomic salt, the renormalized integral equations and the iterative algorithm used in their solution, as well as a specialized Fourier transform have been described in detai1.38v44 We present here the corresponding renormalization for an infinitely dilute model polyelectrolyte. For the present purposes, we utilize a simplified notation in which eq 6 for the rod/ion system is written
h’= c ’ + pc‘* h
Note that the convolution integral is indicated by an asterisk, a sum over small ion species is associated with each p, and the primes indicate rod-ion correlations as opposed to ion-ion correlations. The pair potential between the rod and the ion (eq 4) is divided into us, the short-range core repulsion, and uL,the long-range Coulomb interaction, while 4’ is used to indicate -@uL. We add and subtract 4’ to the second c’in eq A1 and use the definition y’ c‘- 4’ to obtain
* h + p4’ * h
+
(A31
Substituting the quantity h - q q for the last h gives h’ = C’ py’ * h p4’ * ( h - 4) p4‘ * (A4)
+
+
+
where q is the Debye-Huckel screened potential between two ions and is given by
where the Debye
KD
is defined by
KD’
= (4re2@/e)Cp,Zj2 1
(‘46)
We now define
h’ = exp[ql - 1
(A14)
Reorganizing eq A12, using the definitions of Q and cancelling terms, we arrive at
+
h’= exp[q’- @us pc’* h
- p4’ * q ] - 1
T’,
and
(A15)
which explicitly displays the correction terms to the Debye-Huckel approximation. Restating the H N C equations for the rod/ion case in a more detailed notation, we have y’ip= exp[-Pui
