Macromolecules 1993,26, 5414-5422
5414
A Dynamical Monte Carlo Model of Polymer Adsorption Stephen M.King**+ ISIS Science Division, SERC Rutherford Appleton Laboratory] Chilton, Didcot, Oxfordshire OX11 OQX,U.K.
Terence Cosgrove Department of Physical Chemistry, University of Bristol, Cantock’s Close, Bristol BS8 lTS, U.K. Received October 15, 1992; Revised Manuscript Received April 1, 1993. ABSTRACT A novel dynamical Monte Carlo computer model designed to investigate the configurational relaxation of homopolymers at the solid/solution interface is described. Ensembles of polymer chains, each up to 99 segments in length, have been generated by n-step random self-avoidingwalks throughout a simple cubic lattice, the basal face of which represented an impenetrable solid surface. To model adsorption, each chain is generated “in solution”at the center of the lattice and then displaced to the surface and allowed to approach its equilibrium adsorbed conformation through a succession of elementary moves operating on randomly selected segments. Each chain was considered to be in an athermal solvent environment, and the effect of surface coverage was simulated by the inclusion of a periodic boundary constraint. When the enthalpy of adsorption is low, it has been possible to observe chains that have initially adsorbed, relaxed, completely desorbed, and then readsorbed. For isolated nonadsorbed chains, root-mean-square end-bend distance data from the model have provided an estimate of the universal scaling exponent, Y, of 0.691,which is in good agreement with the Flory value of 0.6 for polymers in athermal solvents. By introduction of a probabilityof desorption,the model hae also been used to study the evolution of the structureof the adsorbed layer for polymers undergoing both chemical (irreversibIe)and physical adsorption at the interface. Introduction The adsorption of polymers a t the solid/solution interface is an area of considerable practical and theoretical interest.1 This is largely due to the use of adsorbed polymers as dispersion stabilizers or flocculating agents in many modern industrial processes and products. The conformation of an adsorbed polymer can be described in terms of four parameters. These are the adsorbed amount, I?, or surface coverage, 0; a layer thickness, 6; the average bound fraction of polymer segments, p; and the polymer volume.fraction profile normal to the interface, a(%).In fact, given a complete knowledge of cP(z), it is possible to derive all of the other structural parameters. Experimentally these parameters may be determined by a number of techniques: 6 from photon correlation spectroscopy measurements? electrokinetic methods: or ellipsometry$ p from magnetic resonance5 or FT-IR6 measurements; and @(z) from small-angle neutron scattering @ANSI7or neutron reflectometry.8 Theoretical descriptions of polymer adsorption abound but generally fall into two categories. These differ in how they treat intrachain interactions within the interfacial region and also in how the configurational freedom of the polymer chains is accounted for. Mean-field theories are generally successful in describing adsorption over a wide range of molecular weight, polymer concentration, and solvency. Conversely, scaling theories are analytic and their description of adsorption is in terms of qualitative power l a ~ s . * - ~ * 3 ~ Complementing these experimental techniques and theoretical models are an increasing number of ever more sophisticated computer simulations. The advantage of a computer “experiment” is that it permits one to investigate the effect of many different variables on a given system, t
Formerly of the Department of Physicd Chemistry, University
of Bristol.
Abstract published in Advance ACS Abstracts, September 1,
1993.
either singly or collectively, in a much more systematic way than could reasonably be achieved by practical experimentation. In addition, becausethe simulationmust keep track of the coordinates of every ‘atom” in a polymer “molecule”, the behavior of a chain can be examined in much greater microscopicdetail. Furthermore, it is always possible to calculate analogues of the structural parameters listed above. For reasons largely of computational simplicity, many simulations still utilize the Monte Carlo (MC) method, instead of the more realistic atomistic approaches that are now possible. As its name implies, the MC approach attempts to predict the equilibrium quantities of interest by randomly generating a representative subset of the total number of possible conformations of a chain of a given length. It is a development of this method which we have used in this work. In many simulations of adsorbed polymers it has become customary to “grow” each chain in the ensemble from a point of attachment at the surface and not to let the chains desorb, This is done irrespective of whether or not the simulation is designed to investigate physically adsorbed or grafted chains.l0-l2 In the case of physically adsorbed chains, the reason for adopting this approach is that it reduces the total number of conformations to be sampled. However, in doing so it presupposes that there is no equilibrium between the grafted chains and any polymer in solution. A glance a t a typical polymer adsorption isotherm shows that this condition is, in reality, only met at very low solution concentrations (though the situation is not as clear-cut in the case of polyelectrolytes). Another drawback with this methodology is that because it does not mimic the whole adsorption process only long-time (or what have become known as “static”)ensemble averages are obtained. An alternative procedure is to grow the chains but then move them, and this is the basis of the dynamical Monte Carlo (DMC) method. Although computationally more intense than static MC, DMC is beginning to find favor.
0024-929719312226-5414$04.00/0 0 1993 American Chemical Society
Macromolecules, Vol. 26, No.20, 1993
A Dynamical Monte Carlo Model of Polymer Adsorption 5415
Recent work has focused on polymer adsorption at the liquid/liquid interface,l3polyelectrolyte adsorption at the solid/liquid interface,l4 polymer adsorption on heterogeneous surfaces,15 and chemisorption.49 In a previous paper we have presented bound fraction data from a DMC computer simulation which showed qualitative similaritieswith comparable experimental data from the system poly(methylsiloxane)/alumina/tetrachloromethane.lB In that paper we adopted a novel approach that allowed us to study both the evolution of the structure of the adsorbed polymer layer and the “equilibrium” adsorption. The purpose of this paper is to describe that simulation in detail. Monte Carlo Procedure Our simulation attempts to mimic the configurational relaxation of homopolymers at the solid/solution interface by a two-stage process. First, a polymer chain of up to 99 (n) segments is generated by a n-step random self-avoiding walk (SAW) on a simple cubic lattice. Each segment of the chain occupies one lattice site, and all lattice sites not occupied by the chain are assumed to contain a solvent molecule of equal volume. The first segment is always placed at the center of the lattice, and the remainder of the chain is then generated by the sequential addition of new segments to this initial one. The algorithm checks each new step in the SAW to ensure that it does not contravene the principle of excluded volume. If it does, then that step is rejected and another alternative placement of that step is chosen at random. If no permissible alternative placement can be found within some fixed (but arbitrary) number of attempts, then that chain is rejected and the whole process begins over with the placement of a new first segment. The probability of generating a particular conformation by this process is, however, not the same as the true probability of that configuration occurring on a purely stochastic basis. In turn this introduces a considerable bias into the ensemble which we remove through the system of chain weighting developed by Rosenbluth et al.I7 This ascribes each chain in the ensemble a chain weighting factor, Wn, generated from the product n stepwise transition probabilities, Pi. n
wn= pi t-1
Any ensemble average quantity is then calculated in such a way that it is weighted by Wc
Wc = Wn exp(n,&x,) exp (-xn,/8) and normalized to w
cwc
(2)
N
w=
(3)
c=l
where n d is the number of adsorbed Segments, np is the number of polymemolvent contacts for chain c in an ensemble of N chains, x is the Flory-Huggins polymersolvent interaction parameter ( x = 0 throughout this work), and xais the polymer-surface interaction parameter.l* The quantity 8 is the lattice coordination number. Thus, for example, the ensemble-averagesegment density in some lattice layer z, p ( z ) , is given by C
(4)
where n, is the number of segments of chain c in layer z.
p ( z ) is in fact directly related to the volume fraction of segments, W ) ,through the adsorbed amount, I?. In this way all chains in the ensemble may be counted with equal probability. This method of “choosing” each chain with some Boltzmann-like probability but then weighting it equally, rather than choosing each chain purely at random and weighting it by some Boltzmann-like exponential factor, was first proposed in a general form by Metropolis et al.19 and is now often referred to as Metropolis MC. To minimize the possibility of the program generating the same ensembles of chains on different occasions, the seed for the random number generator was always initialized with a value based on the computer’ssystem clock. Thus, provided that the total number of random numbers requested during a simulation was less than the square root of the period of the pseudo random number sequence used (typically this period is 232 257), the ensemble was unique.
-
The Model The basal face of the lattice (the plane z = 0) is taken to represent an impenetrable solid surface for which the enthalpy of adsorption, xs, can be varied. In this paper xs is always greater than the critical value necessary for adsorption to take place, x,.~O Each newly-generated chain is displaced in the direction normal to the surface (the z-direction) until the segment(8) with the lowest z-coordinate(s) reach the surface. At this point the chain is deemed to be adsorbed. This is the start position for the second stage of the simulation. By generating each chain at the center of the lattice in this way, that chain is free to adopt a “solution-like” conformation. During the running of the simulation an ensemble of N, typically -0(104),chains is accumulated. Each chain is allowed to approach its equilibrium adsorbedconformation through a succession of K elementary moves operating on randomly selected segments. There are a large number of algorithms which could conceivably be used to model the motion of a polymer chain receiving random collisions from solvent molecules. Most use a finite repertoire of local, n-conserving, elementary moves where the term local indicates that the algorithm applies to a short, contiguous, run of chain segments. However, nonlocal, n-conserving, algorithms do exist. One example is the Pivot algorithm.21 In this work we have implemented the set of elementary moves fiist proposed by Verdier et alea These are a trans-gauche end-of-chain move, a gauche-trans end-of-chain move, and a 180’ flip of a pair of gauche bonds known as a “kinkjump”. These moves are illustrated in parts a-c of Figure 1,respectively. We have also included an additional local move, known as the “crankshaft” motion, which has recently been shown as necessary to establish the true dynamics of lattice models of this type.2312~This move involves a 90’ or 180’ flip of two adjacent bonds in gauche positions to the main chain and is depicted in Figure Id. A set of logical rules determines if a randomly selected segment in the chain can move and, if so, what type of move can take place. These rules take into account the position of the segment in the chain, the local lattice occupancy, and the number of nearest-neighbor contacts (i.e., the energy) of both the existing chain configuration and the proposed ~onfiguration.~5 If, during a movementstage, a segment is selected which is already adsorbed, then clearly a decision must be made
6416 King and Cosgrove (8)
T
t
I
0
2*10’
I
I
6*10’ Ncmber of moves. K
4r10’
I
1
8tlO’
lt10‘
Figure 2. Probability of desorption, P d , as a function of the number of elementary moves, K. % = 0.4 (Solid line) Chemisorptivemode. ‘(Dashed line)Physisorptivemode. (Dot-dashed line) Composite mode.
(d)
o--o
Figure 1. Elementarymoves used in the simulation: (a) transgauche end-of-chain move; (b)gauche-trans end-of-chain move; (c) kink-jump move; (d) crankshaft move. about whether to allow that segment to desorb or not. The simulation caters for one of three modes of polymer adsorption. The first of these is a purely “chemisorptive” approach where an adsorbed segment cannot desorb under any circumstances. The second approach is a “physisorptive” one. This allows a segment to desorb only if a probability of desorption, Pd, given by Pd = A + B exp[-(x, - x,JI
(5)
where
A = 11- ~(X,~X,,)IU(X,,-X,) exceeds a random number between zero and 1. In eq 5 , 6 is the usual Kronecker delta and Uis a unit step function with the properties that U(i-j) = 0 if i < j and U(i-j) = 1 if i 2 j . Clearly, for a neutral or unfavorable surface where xsIx,, p d is unity and an adsorbed segment always desorbs. For all other surfaces p d falls exponentially with increasing xs. The final mode of adsorption represents an initially physisorptive interaction between the chain and the surface that becomes more chemisorptive the more the chain relaxes and which we term the composite mode. In this case
+ B exp[-(x,
+ CYCI
(6) where m is the move index (0 5 m 5 K) and C is a constant very much smaller than K and is chosen so as to make p d 0 as m K over the simulation period.26 In this way desorption is made increasingly less likely as m (viz. time) increases. The differences between the three modes of pd
- -
=A
- x,,)(m
adsorption in terms of the variation of Pd with K are shown in Figure 2. The xs and xsc values are the same in each case. The use of an exponential form for p d in eqs 5 and 6 is a purely empirical choice on our part, but a similar idea has been suggested else~here.~7 The simulation incorporates the effect of surface coverage on the adsorbed conformations through a simple periodic boundary constraint.% Every chain in the ensemble is both generated and moved subject to this constraint. Results and Discussion (1) Static Properties. The variations of the ensembleaverage root-mean-square end-bend length, r, and ensemble-average root-mean-square radius of gyration, s, with the chain length, n,for fully isolated, nonadsorbed chains were used to estimate the universal scaling exponent, v, in the expression rasan” where v=[3/(D+2)] and D is the number of dimensions available to the chain.29 For chains with lengths between 7 5 n I63, values of ~0.591 and ~ 0 . 5 4 were 2 obtained, respectively. These estimates are in quite good agreement with the value of 0.6 expected by Flory for homopolymers in a good solvent environment (x < O.5bB Our estimates are also comparable with those of other workers.17v22*30J1 The small discrepancies, most pronounced in the value of v determined from s(n),are attributed to the very small sample size of only N=500 chains in that instance. These results indicated that the simulation behaved as expected and gave confidence that the chain generation routines were satisfactory. The SAW algorithm was further verified by a visual inspection of the lattice coordinates of some chains as output by the program. (2) Effect of Elementary Moves. The next step in the testing procedure was to introduce the elementary moves. Table I details r and s for chains with lengths between 7 5 n - I 63 before and after K-1000 elementary moves. The sample size is again N=500 chains. The data show that r and s are, within the limits of statistical error, unaffected by the incorporation of the elementary moves. The data for K=lOOO give values of ~ 0 . 6 2 2and ~ 0 . 5 6 4respectively. , This is not an unreasonable finding. It could be argued that 1000moves might
Macromolecules, Vol. 26, No. 20, 1993
A Dynamical Monte Carlo Model of Polymer Adsorption 5417
Table I. Root-Mean-Square End-to-End Length, r, and Root-Mean-SquareRadius of Gyration, s, as a Function of the Chain Length, n,and the Number of Elementary Moves,
K*
n
7 15 23 31 39 47 55 63 255
K
r
0 1OOO 0 lo00 0 lo00 0 1OOO 0 1OOO 0 lo00 0 1OOO 0 1OOO 1OOOOO
3.29f 0.56 3.29f 0.64 5.75 f 1.21 5.42 f 1.20 6.91 f 1.57 6.99 f 1.70 7.89 f 2.25 8.26 f 2.11 8.67 f 2.37 9.26 f 2.35 9.55 f 2.55 10.04 f 2.62 10.11 f 2.73 10.64f 2.86 11.09 f 3.04 11.35f 3.07 17.87 f 2.41
9 .
1.40 f 0.13 1.40 f 0.15 2.22 f 0.27 2.17 f 0.26 2.71 f 0.33 2.78 f 0.37 3.11 f 0.42 3.22 f 0.47 3.41 f 0.50 3.58 f 0.53 3.75 f 0.57 3.89f 0.58 3.99 f 0.62 4.12f 0.61 4.34 h 0.70 4.36 f 0.68 7.90 f 1.23
+Is2
5.52 f 1.07 5.52 f 1.23 6.71 f 1.63 6.24 f 1.57 6.50 f 1.68 6.32f 1.75 6.44f 2.03 6.58 f 1.94 6.46 f 2.00 6.69 f 1.97 6.49 f 1.99 6.66 f 2.00 6.42 f 2.00 6.67 f 2.05 6.53 f 2.08 6.78 f 2.12 5.12f 1.05
'Also listed is the constant relating r2 and s2. Quoted errors represent one average deviation.42
not be significant in terms of any perturbation it might inflict on the chain statistics. To test this, simulations with much greater numbers of moves were conducted on ensembles of 99 segment chains. The results are summarized in Table 11. To reduce the computer time necessary, the sample size is only N = 100 chains in each case. While there is a small increase in both r and s with increasing n, once again within the statistical error it would seem that the chain dimensions are unaffected by incorporation of the elementary moves. Table I also shows the constant of proportionality between rZ and s2as determined by the model. The average value is 6.33 f 0.48. In the limit of high molecular weight (equivalent to large n) polymers this constant is expected to tend to a value of 6.00 and so our data shows reasonable agreement. (3) Adsorbed Chains. (3a) Physically Adsorbed. In this section the effect that an impenetrable interface has on the chains and how it affects their conformations is considered. Figure 3 shows the mean initial bound fraction (given the subscript i) as a function of n for isolated chains with lengths between 9 In I99. A solid line has been added to aid the eye. The quantity pi is the instantaneous average bound fraction that exists the moment the newly-generated chain is placed at the interface and is therefore independent of xB and the number of moves, K . The data are highlyrepeatable within the statistical errors, being independent of the choice of random number seed. Each point represents an ensemble of N=20 000 chains. The graph shows that pi decreases with increasing n,a trend that is well described by the empirical power law pi = 2.553r1-~ where y=0.869. This trend may be explained by the fact that, because the chains are isolated, when n is small they are able to adsorb in rather flat conformations which give rise to large bound fractions. However, as n increasesthis is no longer possible because in a random coil the segment density is concentrated around the center of a polymer "sphere" rather than at its periphery. Fewer segments are initially adsorbed as a result and thus pi decreases. It is interesting to note that the data are consistent with the number of initially adsorbed segmentsdoublingwith every 200-fold increase in chain length. How exactly p should scale with n for an isolated chain is unclear. Neville table estimates of y for a cubic lattice
Table 11. Root-Meamsquare End-to-End Length, r, and Root-Mean-Square Radius of Gyration, s, as a Function of the Number of Elementary Moves, K, for 99 Segment Chains. K r S 20 OOO 15.00 f 3.79 5.74 f 1.00 30 OOO 14.80 f 3.59 5.75 f 0.93 40 OOO 15.18 f 3.09 5.82 f 0.86 50 OOO 15.37 f 3.23 5.91 f 0.87 60 OOO 15.97 f 3.15 6.09 f 0.91 70 OOO 15.90 f 3.03 6.08 f 0.88 80 OOO 16.23 f 3.10 6.19f 0.93 16.52 f 2.75 6.35 f 0.85 90 OOO 100 OOO 17.15 f 3.69 6.41 f 1.10 Quoted errors represent 1 average deviation.42
have yielded y=0.89 f 0.02 by n=13 for x s = l but suggest that 7-1 as n--.32 Using some simple geometrical ideas to describe the contact volume between a polymer sphere (the segment density distribution of which is Gaussian about the center) and a planar interface, we have derived a scaling relationship of the form (see the Appendix) p a n-2v Thus, in an athermal solvent y=1.2, decreasing to y=1.0 in a 8 solvent. In our simulation the temporal relaxation or rearrangement of an adsorbed chain is controlled by three parameters: K , xs, and the surface coverage, 8, as implemented by the periodic boundary. For any given value of xs and providing the ensemble is large enough, it is reasonable to expect that the bound fraction will tend toward some limiting value if plotted as a function of (simulation) time or, as in this model, as a function of K . Figure 4 shows an example of one such "convergence test" for ensembles of 49 segment chains adsorbing on an attractive surface (xs=0.82,N=10 OOO chains). It can be seen that even after K=105 moves the mean limiting (or final-hence the subscript f) bound fraction, pf (circles), is still increasing, though the rate of increase could be slowing as the data points for K=105 moves are actually for the slightly higher adsorption energy of x p l . 0 . These required about 20 h on a CRAY supercomputer." The initial bound fractions for each ensemble, pi, are shown as squares. Solid lines have been added as a visual aid. Figure 4 suggests that well in excess of lo6 moves would be required before a (physically adsorbing) 49 segment chain could reach something approximating true thermodynamic equilibrium. Since this work would then have required excessive amounts of computer time, many of the results that follow are from shorter chains where adsorption equilibrium is attained far more rapidly. However, this does not invalidate our results for two reasons. First, it has been demonstratedthat, as equilibrium is gradually approached, 8 only changes by 1-2 3'% per 106 moves at most.15 Second, it must be remembered that the attainment of adsorption equilibrium by real polymers is a very slow process. So slow, in fact, that many experimental studies are probably conducted on nonequilibrium systems. Figure 5 showsthe variation of pfwith xs for an ensemble ofN920 OOO isolated chains of n=9 segments (circles)after K = 1000 moves. It is seen that pf increases with increasing xs as might be expected; the surface becomes more attractive and therefore adsorption becomes more favorable. A dashed line has been used to join the data points as a guide. A similar graph for a cubic lattice has been generated by exact enumeration methods and is depicted in Figure 5 as the continuous line.33 Due to the nature of those calculations, it was possible to extrapolate p(xs)in the limit n -.
-
Macromolecules, Vol. 26, No. 20, 1993
5418 King and Cosgrove
0.~1
‘.O1
0.354
c
3 0.25-
/
t; 0.6
e
S
L
I
/
o.2-
.-.-c
0.15-
=
0.1-
------L50.00.0 0
20
40
60
80
100
I
0.0 0.0
I
0.1
I
I
I
0.2 0.3 0.4
Chain length, n
Figure 3. Mean initial bound fraction, pi, as a function of the chain length, n (N=20 OOO).
I
I
I
0.5 0.6 x. (kT)
I
I
0.7 0.8 0.9
1
1.0
Figure 5. Mean final bound fraction, pt, ae a function of the surface interaction energy, XI: (circles) n=9, (squares) n=49, (triangles) n=99, (continuous line) exact enumeration data for cubic lattice extrapolated to infinite chain length= (N-20 OOO, K=lOOO).
0.51
a 0.41
0.25 Os31
.-$0.3e .u
d
I
L
Q
5 0.2-
d: 0.1 -r
0.0
0
1
2 3 4 LOG (Number of moves, K)
5
6
Figure 4. Mean bound fraction ae a function of the number of elementary moves,K: (squares)initial,pi; (circles)f i i , p r ( ~ 4 9 , N=10 OOO, and x.30.82 except for log (&=5 where xb=l.O).
For a lattice system the value of
xsc is given by36
-.-
I
0.0
0.5
I
I
I
1.0 1.5 2.0 Swfma coverage, 0
I
I
2.5
3.0
Figure 6. Mean final bound fraction, pf, ae a function of the surface coverage, e (n=49, N33000, K=5000, x,=0.82).
xsc= -ln(l-Al)
abscissa is expressed as segments per surface lattice site, 8, such that
The effect of surface coverage on p is shown in Figure 6 for n=49 segment chains adsorbed at xp0.82 after K=5000 moves. The sample size is N = 3 0 0 0 chains. The
8 = n/W2 (8) where W is the width of the lattice (the same in both the x and y directions) in lattice units. Thus for n=49 the chain may be regarded as isolated if 8 < 0.02.9s Reducing the width of the periodic boundary raises the surface coverage. It can be seen that pf falls quite quickly between 0.05 < 8 < 0.2 but then decreases slowly beyond 0-0.68. This latter point approximately coincides with nl(%)2 for this chain length (the actual value would be 8=0.84). Similar observations were made for some lower xs values, but the change in pf was less distinct. Since 2s is the separation at which neighboring chains (in solution) would start to experience one another, the observed behavior can be explained in terms of steric crowding effects. At low 8 the chain is free to relax in such a way that it may maximize ita contacts with the surface and form quite expanded conformations in the plane of the surface. However, as 8 increases, the chain is forced to adopt conformations with more “loops”and “tails” (which extend further away from the surface). With fewer and fewer segments in “trains”, pf is gradually reduced. This behavior has also been seen in MC studies of terminally-attached chains3’ and in
(7) where XI is the fraction of nearest neighbors that could lie in an adjacent lattice layer. Thus for a cubic lattice xsc -0.18, and it is seen that this approximately coincides with the region in whichp(x&+Oforthe exact enumeration data. Although we did not explore the behavior of the simulation in the region xliI xsc,it is clear that our bound fraction data are tending to some small but finite value and not to zero. However, this is also to be expected because xsc is defined in the limit n--.20 To illustrate this point, Figure 5 also shows two additional data points for ~~‘0.22.These are identical data to the circles but for chains of lengths n-49 (squares) and n=99 (triangles). It can be seen that, for this value of xs, as n increases, pf decreases. We do not show any further data points for these longer chains because while lo00 elementary moves might allow a 9 segment chain to reach an equilibrium adsorbed conformation this number of moves would not be sufficient for the longer chains. At very high xSthe exact enumeration data also converge toward a value of p which is close to unity. We have not tested this.
Macromolecules, Vol. 26, No. 20, 1993
A Dynamical Monte Carlo Model of Polymer Adsorption 5419
*lo"
251
r0.5
1
3.0
I
$P 1
20-1
?? 2.5 v a
10.4 2
E
10.15
I
t I
0
1
2
4
3
5 6 z (layers)
7
8
9
10
0
2*10'
4*10'
6*10' "her of moves, K
8*10'
0*05
,
0.0 1*105
Figure 7. Normalized segment density profiles for two different valuesof the surface coverage,8: (circles)low, 8=0.012; (triangles) interpenetrating, 8=2.041 (n=49,N=3000, K=5000,x.=0.82).
Figure 9. Number of adsorbed segments, n h ,and cumulative bound fraction as a function of the number of elementary moves, K, for a chain adsorbing in the physisorptive mode for xpO.8 ( ~ 4 9N=l). ,
2 204
bound fractions after move m in each case, p ( m ) ,as given by
I
0
2*10'
I
4.10'
I
6*10' M e r of moves, K
I
8*10'
I 0.0 1*10'
Figure 8. Number of adsorbed segments, nab,and cumulative bound fraction as a function of the number of elementary moves, K, for a chain adsorbing in the physisorptive mode for xl=O.l 01349,N=l). exact enumeration studies.20~~~39 Figure 6 also compares favorably with predictions from contemporary mean-field theories4 and experimental data.4 Figure 7 shows the normalized ensemble-average segment density profiles for n=49 segment chains at surface coveragesof 8=0.012 and 8=2.041. As the surface coverage increases, it is possible to see a shift in segment density away from the surface and, although less distinct, a small increase in the maximum extent (equivalent to 6) of the profiles. (3b) Chemically Adsorbed. The variation of pf(@has also been investigatedfor the case of chemisorption,where segments are irreversibly adsorbed. It was found that pf was not only smaller than in Figure 6 (for chains of the same length) but also only weakly dependent on 8. This suggests, as might be expected, that the chains are adsorbing in a solution-likeconformationbut then finding themselves unable to rearrange at the surface. (4) Effect of p d . One of the important aspecta of this simulation is its ability to allow segments to desorb, and this is illustrated in Figures 8-10. Each figure shows the number of adsorbed segments, nab, as a function of K for three different (i.e., N = 1) initial conformationsof isolated n=49 segment chains adsorbing under slightly different conditions. Also shown are the corresponding cumulative
In Figure 8 the chain has physically adsorbed. This means that P d is constant and only determined by XI. In this example, xS'0.1, giving Pd'O.9; that is, there is a 90% chance of an adsorbed segment desorbing if selected and if able to move. It is seen that the number of adsorbed segments gradually builds up over the f i t 16 OOO moves but that after 25 OOO moves it completelydesorbs and never readsorbs during the remainder of the simulation (105 moves). In Figure 9 the chain is again purely physically adsorbed but this time the adsorption energy is much higher at x8=0.8. Although this only corresponds to PpO.45, the chain is adsorbed for 1OOO moves, completely desorbs for the next 63000 moves, and then readsorbs with more segments adsorbed than it had originally. Last, in Figure 10the chain has adsorbed with the model running in the composite mode. In this example xs=O.3, thereby allowing Pd to exponentially decrease from 0.74 at K=O to 0.0004 at K=10 OOO. As a result there is an irregular but consistentincrease in the number of adsorbed segments until at the end of the run 44% of them are adsorbed. (5) Segment Density Profiles. XI plays an important role in determining the degree to which a chain will adsorb and it therefore influences the form of p ( z ) . This is illustrated in Figure 11 where p(z), for three ensembles of isolated n=49 segment chains adsorbed at xE=0.22, ~ ~ ~ 0 .and 4 2xa=0.62 , (N=20 OOO chains, K=lOOO moves) in the physisorptive mode are shown. With increasing Xu the maximum in ~ ( 2 diminishes ) ~ and the segment density beyond the first layer drops sharply. The pronounced maximum in the data for xS'0.22 occurs because the system is only 0.04 above xsc. The surface therefore appears almost neutral to the chain. With little to be gained thermodynamically by adsorbing, the chain retains much of its solution conformation. Similar profiles have been obtained experimentally for terminally-attached polymers using SANS.' The surface will actually be even more unfavorable than the preceding discussion impliesbecause, as has already been mentioned, eq 7 is defined in the limit of large n. For shorter chains
Macromolecules, Vol. 26, No. 20, 1993
5420 King and Cosgrove
0’71
25!
-1 0.6
2 20-1
Q
0
I
I
I
2tlO’
4*10’
6+10’
”her
1 0.0
I
&lo’
l*lO‘
of moves, K
Figure 10. Number of adsorbed segments, n h ,and cumulative bound fraction as a function of the number of elementary moves, K, for a chain adsorbing in the composite mode for x,=0.3 (n=49, N=l). 0.71
0.6 v
P
0
1
2
3
4
5
6
7
8
9
10
z (layers)
Figure 11. Normalized segment density profiles for n=49 segment chains physically adsorbing at different interaction energies: (triangles)~ ~ 0 . 2(squares) 2, xn=0.42,(circles)xn=0.62 (N320 OOO, K=1000). the Yeffectivenxscof the surface ri~es.~3 The shapes of the profiles change rapidly as xs increases. Although the maximum extents of all three distributions appear to be the same, it must be remembered that the chain lengths depicted are very short. The important observation to be drawn from these data are that changes in xa only affect the conformations of isolated chains at, or in close proximity to, the surface. The differences between the physisorptive and chemisorptive modes of adsorption suggest two very different segment density profiles for isolated chains. Examples are shown in Figure 12 for ensembles of n=49 segment chains adsorbing at x p O . 2 (minimum N=5000 chains, K=1000 moves). It is seen that in the case of the physisorbed chains ~ ( 2 falls ) ~ rapidly from 0.61 to 0.06 over the first two layers and then decreases more slowly out to about z=6. This is because the chains have been able to rearrange themselvesat the surface,forming trains and loops in the process. In contrast, much less surface rearrangement can take place withthe chemisorbed chains as those segmenta which adsorb do not desorb. Consequently the segment density distribution for this ensemble, ~ ( z ) still ~ , retains much of its initial (bulk solution) Gaussian form. This gives rise to the maximum of p(z),=O.13 at 2=3.
0
2
4
i - i - i - i 6 8 10 12 14 z (toyers)
i - l - 7 16 18 20
Figure 12. Normalized segment density profiles for chains physically and chemically adsorbing: (circles) physisorption, (squares) chemisorption (n=49,”15000+, K-1000, xp0.2).
Conclusions A novel dynamical Monte Carlo computer simulation of homopolymer adsorption at a solid interface from an athermal solution has been presented. Estimates of the universal scaling exponent, v, of 0.591 and 0.542 were made from the root-mean-square end-toend length and root-mean-square radius of gyration data, respectively, for isolated, nonadsorbed chains. The simulation has been used to study the evolution of the structure of an adsorbed polymer layer through the average bound fraction of segments at the surface and the polymer segment density profiles perpendicular to the interface. The enthalpy of adsorption could be varied but was always favorable,and the effect of surface coverage was modeled by the use of a periodic boundary constraint. The trends shown by the data from the simulation are broadly in agreement with published results from both experimental investigationsand contemporary analytical the~ries.~~~~~~~ The simulation is able to model chemical adsorption and, by introducing a probability of desorption, physical adsorption also. In the latter case, and when the enthalpy of adsorption is low, chains can adsorb, desorb, and then readsorb. Besides making the simulation more representative of the real adsorption process, this capability has also enabled qualitative comparisonswith an experimental system to be made.16
Acknowledgment. The authors thank Dr. Iain More (&on Chemicals, U.K.) and Prof. Brian Vincent (Bristol) for their interest in this work. S.M.K.thanks Exxon Chemicals for the award of a postgraduate studentship and Dr.John Webster for several useful discussions. The authors also acknowledge the assistance of Ian Cockerill (Bristol) and Doug Spragg (Exxon Computer Centre, NJ). Appendix In this appendix we derive a simple scaling relationship between the bound fraction of polymer segments at a planar interface, p, and the number of segments in the polymer, n. We start by assuming that the volume occupied by a random-coil polymer in solution can be represented by a sphere whose radius is equal to the radius of gyration of the polymer, s (see Figure 13). So far as the initial adsorption is concerned, we need only consider those segments in a spherical shell between radial distances s-a
Macromolecules, Vol. 26, No. 20, 1993
A Dynamical Monte Carlo Model of Polymer Adsorption 5421 (9) Russel, W. B.; Saville, D. A.; Schowdter, W. R. Colloidal Dispersions; Cambridge Monographs on Mechanics & Applied Mathematics; Cambridge University Preee: Cambridge, U.K., 1989; Chapter 6. 110) Clark. A. T.: Lal. M.: T m i n . M. A,: Richardson, K. A. Discues. Faraday Soc. 1976,'15,i89.' (11) Finch, N. A. Ph.D. Thesis, Universityof Briatol,England, 1988. (12) Cosgrove, T.;Finch, N. A.; Webster, J. R. P. Macromolecules 1990.23,3353. (13) B&s, A. C.; Siemasko, C. P.; Lantman, C. W. J. Chem. Phys. 1991,94,1653. (14) Beltran, S.;Hooper, H. H.; Blanch, H. W.; Prausnitz, J. M. Macromolecules 1991,24,3178. (15) Balazs, A. C.; Gempe, M. C.; Zhou, Z. Macromolecules 1991,24, 4918 and the references cited therein. (16) Coegrove,T.;Prestidge, C. A.; King,S.M.; Vincent, B.Langmuir 1992,8,2206. (17) Rosenbluth, M. N.;Rosenbluth,A. W. J. Chem.Phys. 1966,23, 356. (18) Silberberg, A. Faraday Discussion of the Chemical Society; Royal Society of Chemistry: London, 1978;No. 65. (19) Metropolis, N.;Rosenbluth, A. W.; Rosenbluth, M. N.; Teller, A. H.; Teller, E. J. Chem. Phys. 1963,21,1087. (20) Whittington, S. G.J. Chem. Phys. 1976,63,779. (21) Lal, M. Mol. Phys. 1969,17,57. (22) Verdier, P. H.; Stockmayer, W. H. J.Chem.Phys. 1962,36,227. (23) Stokely, C.; Crabb, C. C.; Kovac, J. Macromolecules 1986,19, 860. (24) The attraction of Verdier's approach lies in both its simplicity and its obvious relation to a real polymer chain subjected to sudden, random impulses. From a practical point of view there are only two constraints on any elementary movement algorithm, (1)that it preserves the connectivity of the chain and (2)that it is ergodic. Essentially, the repertoire of moves in use must allowthe chain the possibility of fully exploring all of the lattice space (computer time permitting) so that one c o n f i a t i o n can be interconverted into any of the other poseible configurations. If this is not possible, then the ensemble w i l l be severely limited. Worse, parts of a chain may gradually become "frozen" until perhaps an entire conformation is "locked". Although not widely appreciated, Heilmann (Mat.-Fys. Medd.-K. Dan. Vidensk.Selsk. 1968,37)and subsequentlyverdier (J.Comput. Phys. 1969,4,204)both realized that Verdier's original set of elementary moves would be nonergodic for three-dimensional SAWS longer than n = 18. This nonergodicity is inherept in these algorithms (see Madras & Sokal, J. Stat. Phys. 1987, 47,573)andcannotberemovedsimplybymakingminorchanges to the logical rules governing movement of the segments. The effect of this nonergodicitywill be to introduce systematic errors into the ensemble; however, it has been demonstrated (see Madras & Sokalor Kremer & Binder, Comput.Phys. Rep. 1988, 7,259)that these systematic errors are much smaller than the statistical errors arising out of the MC procedure. Where there is doubt 88 to the suitability of the typea of moves we have used is in simulations designed to extract dynamical relaxation times or to investigate how these scale with the chain length in the presence of excluded volume (see Verdier & Kranbuehl, Macromolecules 1987,20, 1362 and also Crabb, Kovac et al. Macromolecules 1988, 21, 2230). However, as we make no attempt to do this, we feel justified in our use of these moves. To our knowledge, the incorporation of the "crankshaft" move does not alter the range of ergodicity of Verdier's original algorithm. It has been pointed out to us that at high surface coverage the elementary moves we have ueed are not very effective and might be better replaced by a "snake-like" movement as proposed by Wall & Mandel (J.Chem. Phys. 1975,63,4592). This algorithm is, however, still nonergodic under certain conditions but not to the same degree as Verdier's moves. (25) For example, no move can take place if the selected segment is collinearwith those to either side, whilefewer segmentwlvent contacts imply intrachain interactions which are by definition energetically unfavorable. In this latter case, should the proposed conformation have fewer polymer-solvent contacts, then the move is not made. If the proposed conformation has more polymer-solvent contacts than the existing conformation, then the move only proceeds on the basis of a weighted random choice. This allows for a selection between more than one favorable conformation resulting from alternative types of moves. (26) In our previous paper we chose to express eq 6 with a simulated rate constant, k,, and took xr into x,; hence, p d = A + B e x p [-(x. + k,m)l. The two forms of the expression are identical. (27) Balazs, A. C.; Gempe, M. C.; Brady, J. E. J. Chem. Phys. 1990, 92,2036. .
R
*-
Contact Sector
/
Figure 13. Cross-sectional representation of a random-coil polymer of radius of gyration s and segment length a adsorbed at a solid surface. See the Appendix.
and s, where a is the length of a segment. The volume of such a shell, V, is then given by
v = (4/3)?r(3s2a - 3sa2+ a3) If we further assume that there is a sphericallysymmetric distribution of segments about the center of mass of the polymer, then the segment density, p(R), at some radial distance, R, is given b y the well-known expression The number of segments at distance R from the center, NR,is simplyproportional to the product Jp(R)ds V, where the integral is taken between s-a and s. However, since the segment density in this region is small and s>>a, we use the approximation that N R - ~ ( R V. ) When R=s NR C,(n/s3)(3s2a - 3sa2 + a3) where C1 is a constant. Only a fraction, F, of the spherical shell will actually be in contact with the interface. For a spherical sector of thickness a N
F = 2?rsa/4?rs2= a/2s and the number of segments within it, Na, is j u s t the product FNR. Utilizing the scaling relationship sa nuand neglecting terms of O(a3)and higher sinces>>a,it m a y be shown that N , C2n1-2u where C2 is a constant. It then follows that N
p = N,/n
-
C2n-2u
References and Notes (1) Coegrove, T.; Vincent, B.; Cohen-Stuart, M. A. Adu. Colloid Interface Sci. 1986,24,143. (2) Killman, E.;Maier, H.; Baker, J. A. Colloids Surf. 1988,31,51. (3) Cohen-Stuart, M.A.; Tamai, H. Macromolecules 1988,21,1863. ( 4 ) Takahashi, A.;Kawaguchi, M.; Kayashi, K.; Kato, T. ACS Symp. Ser. 1984,240,39. ( 5 ) Cosgrove, T.; Vincent, B.; Cohen-Stuart, M. A.; Barnett, K. G.; Siseons, D.S.Macromolecules 1981,14,1018. 16) Coeer0ve.T.: Preatidee. C. A.: Vincent. B. J. Chem. Soc.. Faraday !I'rgns.1990,86,13f7.' Cosgrove,T.; Crowley,T. L.; Vincent, B.; Barnett, K. G.; Tadros, Th. F. Faraday Symposia of the Chemical Society, Royal Society of Chemistry: London, 1982;No. 16. Rennie, A. R.; Crawford, R. J.; Lee, E. M.; Thomas, R. K.; Crowley, T. L.; Roberta, S.; Qureshi, M. S.; Richards, R. W. Macromolecules 1989,22,3466.
I
6422 King and Cosgrove (28) Clark, k T.; Lal,M.Polymer 1975, 16, 310. (29) de Gennea, P.-G. Scaling Concepts in Polymer Physics; Cornell University Press: Itacha, NY,1979. (30) Kolinski, A.; Skolnick, J.; Yaris, R. J. Chem. Phys. 1987, 86, 7164. (31) Milik, M.; Orszagh, A. Polymer 1990,31, 506.
(32) Coegrove, T. Unpublished work. (33) Cosgrove, T. Macromolecules 1982,15, 1290. (34) The remainder of this work waa conductad on DEC VAX 6210 and GOULD PowerNode (both about 25 times slower) and the author’s DEC MicroVAX 2000 (about 40 times slower) computers. (35) Scheutjens, J. M. H.M.; Fleer, G. J. J. Phys. Chem. 1979,83, 1619. (36) Other workers have expressed the surface coverage in terms of a dimensionleesreduced area per chain (see ref 37, for example),
Macromolecules, Vol. 26, No. 20, 1993 A. This is related to 0 through6=1/(LA) where L is the contour length of the chain (equivalent to n in this model). (37) Clark, A. T.; Lal, M.J. Chem. Soc., Faraday Trans. 1978, 74, 1867. (38) Lax,M. J. Chem. Phys. 1974,60,2627. (39) Torrie, G. M.; Middlemias, K.M.; Bly, S. H. P.; Whittington, S. G. J. Chem. Phys. 1976,65,1867. (40) Cosgrove, T.; Heath, T.; Van Lent, B.; Leermakers, F.; Scheutjens, J. Macromolecules 1987,20, 1692. (41) Cwgrove, T.; Cohen-Stuart, M. A.; van der Beek, G. Langmuir 1991, 7, 327. (42) Press, W. H.; Flannery, B. P.; Teukolsky, S. A.; Vetterliig, W.
T. Numerical Recipies in C Cambridge - Universitv Press:
Cambridge, U.K., 1988. (43) Konstadmidie.. K.:. Praaer, - . S.:. Tirrell, M. J . Chem. Phvs. 1992. 97, 7777.
