DMRG studies of the effect of constraint release on the viscosity of polymer melts

6 Pages·0.16 MB·English

by  Matthias Paessens

DMRG studies of the effect of constraint release on the viscosity of polymer melts
Matthias Paeßens∗ and Gunter M. Schu¨tz
Institut fu¨r Festk¨orperforschung, Forschungszentrum Ju¨lich - 52425 Ju¨lich, Germany
(Dated: February 1, 2008)
2 Thescaling oftheviscosity ofpolymermeltsisinvestigatedwith regard tothemolecular weight.
0 We present a generalization of the Rubinstein–Duke model, which takes constraint releases into
0 accountandcalculatetheeffectsontheviscosity bytheuseoftheDensityMatrixRenormalization
2 Group(DMRG) algorithm. Usinginputfrom Rousetheorytherates for theconstraint releases are

n determined in a self consistent way. We conclude that shape fluctuations of the tube caused by
a constraint release are not a likely candidate for improving Doi’s crossover theory for the scaling of
J thepolymer viscosity.

5
2
I. INTRODUCTION along which the displacement of the polymer chain as

] a whole can be monitored. This makes it possible to
h
calculate the diffusion coefficient and the viscosity from
c The viscosity of polymer melts has been investigated
e intensively.1 The models for describing this behavior, the model,without resortingtoindependent hypotheses.
m Rubinstein assumes that the constraints of the other
however, are not yet satisfying. Experiments show that
- the viscosity η scales like M3.3±0.1, where M is the polymers divide space into cells which form a d–
t dimensional regular cubic lattice. The polymer occupies
a molecular weight. This behavior is valid for several
t decades of the molecular weight. An early approach a series of adjacent cells, the “primitive path”. It is not
s
. was the reptation model by de Gennes2 which yields possibleforthe polymerto traversethe edgesofthecells
at η M3 in the limit of infinitely long polymers. For (in two dimensions: the lattice points) so that only the
m sho∝rt polymers Doi3 calculated the viscosity by taking ends of the polymer can enter new cells. The polymer is
divided into segments whose length is of the orderof the
- tube length fluctuations into accountandfound a region
d of M where the scaling is of the correct size – but this lattice constant, the number of segments is proportional
n to the length of the polymer or the molecular weight.
region appears to be much too small. A discrete model
o The orientation of the lattice is introduced by Duke in a
for reptation which includes tube length fluctuations is
c
[ the Rubinstein–Duke (RD) model.4,5 Recently it was way that the electric field is diagonal to the lattice, i.e.
shown6 that this model does not only provide a region in three dimensions in the (111)–direction.
1 where a scaling of M3.3±0.1 can be found, as already A segment – called a “repton” – is allowed to jump into
v shown by Rubinstein,4 but also that this model shows an adjacent cell according to the following rules:
3
6 a crossover to the reptation exponent 3. But again the
4 mass regionwhere the correctexponent is validdoes not
1. The reptons in the bulk are only allowed to jump
1 exceed one order of magnitude. In this paper we extend
along the primitive path.
0 the RD model in order to investigate whether constraint
2
release (CR) broadens that region. While good models 2. No cellinthe interiorofthe primitivepathmaybe
0
/ for the non–linearregimeof viscositycanbe constructed left empty.
t byusing CR,7,8 therewasonlylittle successinthe linear
a
m regime.9 3. The ends move freely insofar rule 2 is respected. If
an end repton occupies the cell alone, it can only
-
d retract in the cell of the adjacent repton. If the
n adjacent repton is in the same cell, the end repton
o may enter any of the 2d surrounding cells. Rep-
II. THEORY
c
tons inthe bulk jump withthe sameprobabilityas
:
v reptons at the end into occupied cells.
A. RD model
i
X
Rule 1 ensures that the polymer does not traverse the
r In order to investigate the role of tube length fluc- edges of the cells. Rule 2 is motivated by the fact that
a
tuations in reptation theory Rubinstein introduced a the segments are of the size of the lattice constant. Fi-
discrete “repton” model which allows the description nally rule 3 reflects the fact that there are more free ad-
of three–dimensional reptation dynamics by a one– jacent cells for an end repton than occupied ones. If one
dimensional lattice gas.4 Duke generalized the model considers a field F = 0 the ratio of the probabilities to
6
to the case that an external electric field acts on a jump in respectively against the direction of the field is
chargedpolymer5 which allowsa gooddescriptionofthe proportional to the Boltzmann weight. In the following
diffusion constant in gel electrophorese experiments.10 wewillonlyfocusonthecasewithoutfield,weonlyneed
However,also in the absence of a field the generalization the reference axis.
by Duke is useful in so far as it provides a reference axis As the shape of the primitive path is not affected by
2
F 6. At the ends of the chain A–particles may be cre-
ated: ∅ A.
→
1 2
7. At the ends of the chain B–particles may be cre-
3
4 ated: ∅ B.
5 0 B B 0 B 0 B A →
6
9 The processes 1, 2, 4 and 5 take place with the same
7
8 probability; the processes 6 and 7 are d–times more
probable.
In this way a stochastic interacting particle system
a) b)
on a one dimensional chain has been defined. As the
transitions are independent of the previous history this
process is Markovian.
FIG. 1: a) The repton model in two dimensions: the circles
represent the reptons; the primitive path is marked by the
bold lines. b) Projection onto onedimension.
B. Quantum Hamiltonian
the movement of the polymer and only the ends are cre-
ated or annihilated this model can be mapped onto one A convenient way to describe the process mathemat-
dimension. To this end a particle is assigned to each ically is the quantum Hamiltonian formalism which we
bond between reptons; starting on one end, a bond into will present here shortly, for details see Ref. 11.
the direction of the field is identified with an A–particle, The probability to be in state η at time t is labeled by
| i
against the field with a B–particle and finally a bond P (t). These probabilities of the individual states can
η
without change of potential, i.e. a bond between two be combined to a vector: P(t) = P (t) η . Due to
η
reptons in one cell, by a vacancy ∅ (Fig. 1).17 the conservation of probab|ilityithe ePntries o|f tihe vector
These particles are residing on a chain whose number of P(t) sum up to 1 at any time t. With this definition
| i
sites N is the number of bonds, thus the number of rep- the master equation can be written as:
tons minusone. As foreachsite three statesarepossible
d
this model can be identified with a “quantum” spin–one P(t) = H P(t) , (1)
dt| i − | i
chain.11,12 We define
with the stochastic generator
1 0 0
A 0,∅ 1,B 0. H =− wη′→η|ηihη′|+ wη→η′|ηihη|.
≡ 0 ≡ 0 ≡ 1 Xη ηX′6=η Xη ηX′6=η
      (2)
We label this set of states with X = {A,∅,B}. A state Here wη′→η is the transition probability from state |η′i
of a chain of length N can be considered as an element to η . In other words the off diagonal elements of the
η of the tensor base X = X⊗N; η is constructed by ma|triix H are the negative transition rates between the
| i | i
thetensorproductofthethreecomponentvectorsforthe respective states and the diagonal elements are the sum
individual sites. of the rates leading away from the respective state.
The dynamics of this one dimensional model is: Thecreationoperatorsa† andb† aredefinedbya†∅=A
and b†∅=B, acting on any other state yields zero. The
1. A–particlesmayexchangewithvacancies∅: A∅⇌ annihilation operators a and b are defined by aA = ∅
∅A. andbB =∅,againactingonanyotherstateyieldsagain
2. B–particles may exchange with vacancies ∅: zneBro=. Fbi†nbaallnydwne∅de=fin1e thneAnumnbBe.r oBpyertahteosresdneAfin=itaio†nas,
B∅⇌∅B. the stochastic generato−r H of−the RD model reads:12
3. A–particles may not exchange with B–particles. N−1
H =b (d)+b (d)+ u (3)
1 N n
By these rules it is guaranteed that the shape of the X
n=1
primitive path is conserved. The next rules define the
with
boundary dynamics:
b (d) = d n∅ a† +n∅ b† +nA a +nB b
4. At the ends of the chain A–particles may be anni- n n − n n − n n − n n − n
hilated: un = nAn(cid:2)n∅n+1−ana†n+1+n(cid:3)Bnn∅n+1−bnb†n+1
A ∅. + n∅nA a†a +n∅nB b†b .
→ n n+1− n n+1 n n+1− n n+1
5. At the ends of the chain B–particles may be anni- Again, d labels the lattice dimension. We have chosen
hilated: the time scale such that the hopping rate in the bulk
B ∅. equals unity.
→
3
C. Calculation of the viscosity
The viscosity is proportional to the longest relaxation
time of the stochastic generator H which is the inverse
energy gap.3 Due to the conservation of probability the
groundstate ofa stochasticgeneratorhas the eigenvalue
zero, so that we only need to calculate the first eigen-
value.
A very efficient algorithm to calculate the lowest excita-
tionsofquantumspinchainsisthedensitymatrixrenor-
malization group (DMRG) algorithm.13,14 For the stan-
dardRDmodelthe usefulnessofthisalgorithmhasbeen
demonstrated,6 so that we employ this algorithmas well
for the modifications of the RD model, to be introduced
in what follows. FIG. 2: A constraining polymer vacates for a short time its
The Hilbert space of the spin chain grows exponentially site and thusmakes a constraint release event possible.
with the number sites N; the number of states is 3N,
so that for N = 50 — the maximum number of sites F
considered in our calculations — the number of states
is of the order 1023. By this it is reasonable to project
all quantities onto a subspace consisting of the most im-
portant states. The difficulty is to find out which states
are the “most important”. The DMRG algorithm is a
method which provides a choice of states which is opti-
malintermsofamaximumoftheprobabilitywithwhich
the states contribute to the target state.13 This proba-
...AABA... ...ABAA...
bility is gained by diagonalizing the density matrix. As
the diagonalization of the density matrix for the whole
system would be as laborious as the direct diagonaliza-
FIG. 3: Constraint release by AB–permutation in the RD
tion of the hamiltonian the system is build up stepwise.
model.
Starting with a small system (e.g. 2 states) one adds
iteratively two states until the system has reached the
searched size. In each step the Hilbert space is reduced (CR) event.
to the subspace of the most probable states, so that the A straightforward implementation of this mechanism is
systemsizeincreaseswhilethedimensionofthematrices shown in Fig. 3. The turning of bend of the primitive
remain constant. path equals to the permutation of an AB–pair. At the
ends ofthe pathA–particlescanbe transformedintoB–
particlesandviceversawhichcorrespondstoaCRevent
D. Constraint release at the ends.
ToconstructtheextendedHamiltonianwhichtakesthese
The RD model as well as standard reptation theory is mechanismsintoaccountweintroduceanoperatorwhich
based on the assumption that the polymer moves in a transforms B–particles into A–particles:
fixednetworkformedbythesurroundingpolymers. This c=a†b; c† =b†a, (4)
assumption is at best justifiable for a single polymer im-
mersed in a gel. But for polymer melts, for which the the adjoint operator effects the reversed process. The
viscosity is measured, it should be taken into account corresponding diagonal elements are build up by nB re-
that the surrounding polymers reptate themselves. In spectively nA. Now the new Hamiltoniancan be written
terms of the Rubinstein model this means that the lat- as:
tice itself is subjected to fluctuations. N−1
We consider the following model of lattice fluctuations: HCR =H +g (α)+g (α)+ v (α)
1 N n
Imagine that a constraining polymer moves so far that X
n=1
the constraintfor the investigatedpolymer is releasedso g (α)=α nB c +nA c†
that it can move freely in this region. After a short time n n − n n − n
(cid:2) (cid:3)
theconstrainingpolymerreturnsoranotherpolymerhas v (α)=α nBnA c c† +nAnB c†c ,
taken its place so that the free movement in this region n h n n+1− n n+1 n n+1− n n+1i
(5)
isagainprevented. Thelatticehasregaineditsoriginally
structurebuttheprimitivepathmayhavechangedinthe with the Hamiltonian H from Eq. (3), and α labels the
bulk (Fig. 2). This iswhatwecalla “constraintrelease” rate for the CR process.
4
III. CALCULATIONS AND RESULTS
In order to investigate the dependence of the viscosity 10000
on the CR rate α we performed DMRG calculations for
therates1/N,1/N2 and1/N3. Thechoiceoftheserates
isbasedontwoestimates. Thefirstisverysimple: Acon-
straintis releasedwhen anend repton retractsback into
η
the tube. The probability that the constraint is exerted
by an end repton is 2/N, provided that the probability
for the presence of a segmentis equally distributed. The
without CR
hopping ofa singlereptonis a processofrateone andso Rate ~ 1/N3
Rate ~ 1/N2
we find that the rate for a CR is of the order 1/N. For Rate ~ 1/N
the second estimate we assume that a CR is causedby a
tube renewal. Later in this paper we will show the cal- 1
10 100
culation of this estimate within the Rouse model which N
leads to the result that the rate α is proportional to the
inverse relaxation time and scales therefore with 1/N3.
FIG. 4: The viscosity (η ∝E−1) for several CR rates of AB
The viscosity which was calculated with these rates for N
permutation. Thesolidlineshowstheresultforthestandard
the mechanism of AB exchange is plotted in Fig. 4. It
RD model without CR.
can bee seen that the viscosity decreases with increasing
CRratewhichisreasonablebecausetherelaxationofthe
tube is acceleratedby this process. For the investigation
of the scaling it is more useful to plot the local slope or
the effective exponent
lnτ lnτ
z = N+1− N−1 , (6) 3.4 without CR
N ln(N +1) ln(N 1) Rate ~ 1/N3
− − Rate ~ 1/N2
Rate ~ 1/N
against1/√N asintroducedinRef. 6. The choiceofthe
abscissa is motivated by the formula of Doi3 which pre- zN 3.2
dictsacorrectiontotheN3–scalingintheorderof1/√N
bytakingtubelengthfluctuationsintoaccount. Itshould
be mentioned that in this paper N labels the number of
3
segmentswhichisproportionaltothe lengthofthe poly-
merwhileinRef. 6itlabelsthenumberofreptonswhich
isthe numberofsegmentsplusone. Surprisinglythe lat-
2.8
terinterpretationshowsabetteragreementwiththeDoi 0 0.1 0.2 0.3 0.4 0.5
1/N1/2
formula while the former shows a better agreement with
the experiments since the range where z 3.3 0.1 is
N
≈ ±
much broader. FIG.5: TheeffectiveexponentzN forseveralCRratesofAB
The influence on the effective exponent zN is plotted in permutation. Thesolidlineshowstheresultforthestandard
Fig. 5. The influence is non–monotonic in α: While RD model without CR.
the rate 1/N causes an obvious shift down in compari-
son with the data without CR, the data of rate 1/N2 is
located above the curve without CR. The data of rate
1/N3 is also located above the curve without CR but at that the first excitation is caused by an other relaxation
a smaller distance. mechanism and the second excitation yields the desired
We remark that from a theoretical point of view also time, which scales with the expected exponent whereby
the rate α = 1 is interesting. If A– and B–particles ex- the transition to Rouse dynamics is verified.
changewiththesamerateasparticleswithvacancies,the Thebounds1/N >α>1/N3areratherweakandaself–
particle can diffuse freely, as they feel no longer restric- consistentapproachisneededfordeterminingtheappro-
tions. This means that the reptation model transforms priate scaling of the rate α. As a CR event is caused
into the Rouse model as the tube can change its form by the tube renewal of a polymer, this process is not in-
freely. This is why rate α = 1 represents a possibility dependent of the relaxation time. The relaxation time
to verify the model: The relaxation time should scale in itself is affected by the CR rate. So the rate is physical,
the limit N as N2. But the DMRG calculation iftherelaxationtime,whichiscalculatedusingthisrate,
→ ∞
ofthe effectiveexponentwiththis rateyieldsadiverging yields the same CR rate.
curve. Ontheotherhand,ifthesecondexcitationiscon- To calculate the self consistent rate α we proceed as
SC
sidered, one can see that it scales as N−2. This suggests follows: ThedependenceoftheCRrateontherelaxation
5
1e+05
τ(α)
14000 Ν
τ(α)
12000 10000
10000
1000
8000
τ η
6000 100
4000
10 without C.R.
2000 S.C.
DMRG with α
SC
00 0.0005 α 0.001 0.0015 1 10 N 100
FIG. 6: Relaxation time versus CR rate α; αSC is given by FIG. 7: The new DMRG calculation of viscosity – a good
theintersectionpointofthecurvesτN andτ (hereexemplified coincidence with the data of the intersection points can be
for N =30). seen.
time α(τ) is estimatedby ananalyticalcalculationusing is the first passage time density, i.e. the probability per
Rouse dynamics. The DMRG calculation with the rates time that the segment s is reached by one end just at
mentioned above yields directly the dependence of the time t. Now we can specify the mean first passage time
relaxation time on the CR rate τ (α) for the respective
N ∞
chainlengths. Theinversefunctionofα(τ),τ(α),iscom- µ¯(s)= dttf(s,t) (10)
Z
pared with τN(α): The intersection of the curves yields 0
the self consistent rate α in first approximation.
SC whichindicateshowlongittakesinaverageuntilsegment
First we calculate the dependence of the CR rate on the
s is reached by an end. Finally we average over all s
relaxation time. The polymer whose relaxation time we
provided that each segment builds up an entanglement
want to calculate finally is hindered in his free move-
with the same probability
ment by an other polymer. The position which exerts
the constraint may have the distance s from one end of 1 N
the constraining polymer. So, the constraint is released µ= dsµ¯(s). (11)
N Z
when the polymer moves either the distance s into the 0
one direction or the distance N s into the other direc-
− This is the desired time: In averagethe release of a con-
tion, where N labels the length of the polymer. In Ref.
straint will take the time µ, the CR process takes place
3 anexpressionis presentedwhich providesthe distribu- with rate α = µ−1. The calculation of the integrals is
tion ψ(s,t) of the probability that the tube segment s
elementary mathematics and we only show the result:
was at time t not yet reached by the ends:
µ=0.822τ α=1.22τ−1. (12)
ψ(s,t)= 4 sin pπs e−p2t/τ, (7) ⇒
X pπ (cid:16) N (cid:17) This calculation is based on a continuous description
p,odd
of the polymer. However, we are considering a lattice
with the relaxation time τ. The calculation of this for- model of reptation so that boundary effects might
mula bases on the assumption that the polymer diffuses not be taken into account correctly by the continuum
inbetweenthetubeduetoRousedynamics. Thecomple- expression (11). Therefore we performed the above
mentary probability distribution φ(s,t), which indicates calculation as well for the discrete case. As the result
the probability that at time t the segment s is already for the discrete case converges quickly to the continuous
reached by one of the ends is then case but consists of only numerical evaluable series, we
restrictedourselvestothe continuouscasehere. Thereis
φ(s,t)=1 ψ(s,t). (8) anerroroflessthan2percentforaslittle as10reptons.
−
The derivate with respect to time Fig. 6 shows the curve τ(α) = 1.22α−1 and the re-
laxation times for the CR rates α = 1/N3 and 1/N2
∂φ exemplified for N = 30, the data point for α = 1/N
f(s,t)= (9)
∂t is out of the plotted range. In order to investigate the
6
4.2 new DMRG run was done with a fit of the determined
rates; as can be seen in Fig. 7 the resulting relaxation
4 times coincide well with the ones determined by the self
without C.R. consistence condition, so that the linear interpolation
S.C.
DMRG with α represents a sufficiently good approximation.
3.8 SC
InFig. 8the effectiveexponentisshownfortheselfcon-
sistent relaxation times and for the DMRG calculation.
zN 3.6 The exponent is shifted to higher values in the region of
small N, but a broadening of the region where z 3.3
N
3.4 could not be observed. ≈
3.2
The Rubinstein model is not microscopic and hence
3 there is some freedom in implementing microscopic pro-
0 0.1 0.2 0.3 0.4 0.5
1/N1/2 cess such as CR. A different mechanism10 leads to the
creation and annihilation of particles in the bulk. We
performedsimilaranalysisforthis mechanismwithqual-
FIG. 8: The effective exponent zN of the viscosity with self– itatively similar results.16
consistent CR.
We conclude that Rouse–based calculations do not lead
to a proper description of crossover behavior in the vis-
cosity of polymer melts when Rouse theory is used as
dependence of the relaxation time on the CR rate we a self consistency input in the mesoscopic and generally
interpolated the data points linearly by fitting them quitesuccessfulRubinsteinmodelforreptation. Wecan-
with splines. The data points are is well approximated not rule out that a fully self–consistent implementation
by this interpolationinthe regionwhere the intersection of CR (without Rouse assumption) leads to a crossover
point is expected. In this way we determined the self regime closer to empirical evidence, but the broad range
consistent relaxation times and rates for the lengths of CR rates studied here suggests that bulk shape fluc-
N = 8...50. As the lengths N = 4 and N = 6 do tuations of the tube which result from constraint release
not lead to reasonable intersection points we did not do not significantly broaden the crossover range of the
take them into account. To verify the self consistence a viscosity.
∗ Electronic address: [email protected] 10 G.M. Schu¨tz, Europhys.Lett 48, 623 (1999).
1 J.D. Ferry, Viscoelastic Properties of Polymers (Wiley, 11 G.M. Schu¨tz, in Phase Transitions and Critical Phenom-
New York,1980). ena, editedbyC. DombandJ.Lebowitz(Academic,Lon-
2 P.G. de Gennes, Scaling Concepts in Polymer Physics don, 2000), Vol. 19.
(Cornell University Press, London,1979). 12 G.T. Barkema and G.M. Schu¨tz, Europhys. Lett. 35, 139
3 M.DoiandS.F.Edwards,TheTheoryofPolymerDynam- (1996).
ics (Clarendon Press, Oxford, 1986). 13 S.R.White, Phys. Rev.Lett. 69, 2863 (1992).
4 M. Rubinstein,Phys.Rev.Lett. 59, 1946 (1987). 14 R.M.NoackandS.R.White,inDensity–MatrixRenormal-
5 T.A.J. Duke,Phys. Rev.Lett. 62, 2877 (1989). ization, A New Numerical Method in Physics, edited byI.
6 E. Carlon, A. Drzewinski, and J.M.J. van Leeuwen, Phys. Peschel and M. Kaulke(Springer, Berlin, 1999).
Rev.E. 63, 010801(R) (2001). 15 M.HenkelandG.Schu¨tz,J.Phys.A:Math.Gen.21,2617
7 C.C. Hua, J.D. Schieber, and D.C. Venerus,J. Rheol. 43, (1988).
701 (1999). 16 M. Paeßens, Diplom thesis, RWTH–Aachen (2001).
8 A.E. Likhtman, S.T. Milner, and T.C.B. McLeish, Phys. 17 In the Rubinstein model only two cases are distinguished;
Rev.Lett. 85, 4550 (2000). a bond between two cells is identified by a particle and a
9 N.A. Rotstein, S. Prager, T.P. Lodge, and M. Tirrell, bond within a cell is identified bya vacancy.
Theor. Chim. Acta 82, 383 (1992).