353
Journal of Non-Newtonian Fluid Mechanics, 2 (1977) 0 Elsevier Scientific Publishing Company, Amsterdam
A NEW CONSTITUTIVE THEORY
NHAN PHAN Department (Australia) (Received
THIEN
and ROGER
of Mechanical November
EQUATION
353-365 - Printed
in The Netherlands
DERIVED FROM NETWORK
I. TANNER
Engineering,
University
of Sydney,
Sydney,
N.S. W. 2006
23,1976)
Summary
A constitutive equation is derived from a Lodge-Yamamoto type of network theory for polymeric fluids. The network junctions are not assumed to move strictly as points of the continuum but allowed a certain “effective slip”. The rates of creation and destruction of junctions are assumed to depend on the instantaneous elastic energy of the network, or equivalently, the average extension of the network strand, in a simple manner. Agreement between model predictions and the I.U.P.A.C. data on L.D.P.E. is good.
1. Introduction
Several constitutive equations derived from the network picture of polymer solutions have been proposed [ 1,2]. These models are not completely satisfying in the sense that either the final stress tensor is expressed in a configuration integral involving microscopic variables or material responses are not modelled accurately due to inadequate assumptions about the rates of creation and destruction of network junctions. In two papers by Wiegel [3] and Wiegel and De Bats [4] which seem to be overlooked in the network literature, the problem was approached along the line of Boltzmann’s kinetic theory of gases. The stress tensor was shown to assume a Boltzmann integral form where the kernel was specified in terms of the creation and destruction rates. In this paper we carry the Wiegel formulation a step further by not assuming an affine transformation relating strand stretching and gross motion of the fluid. Specific forms of creation and destruction rates are assumed based on the idea that the instantaneous elastic energy of the network strand is impor-
354
tant. A constitutive equation for the stress tensor is derived. Agreement Meissner’s [ 51 data on L.D.P.E. is good to excellent. 2. Problem
with
formulation
Throughout this paper we will be using the network terminology proposed by Lodge [ 21. A typical example of the kind of network we shall be concemed with is shown in Fig. 1. Each network strand is represented by vector p between its endpoints. We will neglect the distribution of N, the number of subunits between two junctions, for the present. The effect of N will be allowed for later on by assuming multiple relaxation times. We will be considering only homogeneous flows of incompressible isothermal polymeric fluids. The velocity gradient tensor L thus satisfies trL=O Consider the motion We have
(1) of a junction
having position
r(t + At) = r(t) + LrAt + aAt, where OLis certain “slip vector”.
r(t) over a short period At. (2)
Hence
P(t + At) = r2(t + At) - rl(t + At), p(t + At) = p(t) + LpAt + (a(~) - a(rl))At.
Fig. 1. Typical network of polymer solutions.
355
To first order in Ipi we can write b(t) = Lp + g
p .
(3)
For acr/ar = 0 we arrive at the customary assumed “affine transformation” [ 21, which assumes that the network and fluid strains are the same about a point. Now it is reasonable to assume that the “slip tensor”
(4)
where pj are functions of the appropriate invariants of D, p combinations. Equation (4) is seen to be identical to Ericksen’s anisotropic fluid theory [7] where the fluid is assumed to be composed of single-directive segments. To simplify the mathematics we will assume the following non-affine transformation, which is a special case of (4), B= Lp-up--tDp,
(5)
where u, t are constants. Gordon and Schowalter [8] have used an identical expression in their linear bead-spring model. Now let f be the probability distribution of junctions such that f d3 p d3r represents the number of junctions in the statistical volume d3 p d3 r. The total derivative of f is af af at+v.-ar+I,
.$=(
ra t e o f creation tions.
- rate of destruction)
of network
junc-
Furthermore we assume that the rates of creation and destruction of junctions are represented by g(p) and h(p)f respectively [ 1,2,4]. The conservation equation for f becomes
af,, .,+(Lp--op--Dp)+=g--hf. af
(6)
at
(6) has been derived by Wiegel [ 31 using the affine transformation assumption in which u = 5 = 0. Assuming a Gaussian network, the stress tensor is given by the convolution integral [l-3]
T-3)20 -iiF
ppfd’p=;$ “UY
[PPI 9
(7)
where 0 is the absolute temperature, k is the Boltzmann constant, (IVCZ~)~‘~ is the equilibrium extension of the network strand and [ ] denotes averaging over the configuration space (p-space).
356
The problem reduces to solving equation (6) for the distribution function f. The stress field is then found by the convolution integral (7). 3. Constitutive equation We can avoid solving eqn. (6) for f by multiplying it by pp and averaging both sides over the configuration space. By assuming that lim p”f= 0 P-‘” it can be shown that $ [PPI -Z =
jjjaw
[PPI - [PPI LT - 8ZJ]PPl
d3p- jjjhwf
+ [PPIZJ) + 5 c
[PPI (8)
d3p.
Now the elastic (Helmholtz) energy of the network strand is simply proportional to the average extension of the strand [ 21, [p’]. If we take this to be the important parameter and assume
(9) h = WP21),
(10)
then equation (8) thus becomes dT --dt
.C’T-T.f?T
+HT=NkeGZ,
(11)
where Z is the identity tensor, P = L - [D is the effective velocity gradient and H( [p2]) is simply h( [p2]) + 5 u. It is not absolutely necessary to proceed to the result (11) in this way. For example, dropping the averaging inside the G function in eqn. (9) will still yield an isotropic tensor. In the integral for H we can also see that the proposed result is true in the extreme cases when h is not much affected by p (“weak” flow) and also when a very “strong” flow exists [9] so that the distribution function f is practically a delta function localized in p-space; thus we expect the approximation adopted to be useful over a wide range of flow rates. By writing eqn. (11) in convected form [lo] we can solve for the stress tensor explicity : T-
Nke j G(t’) exp{-aa
_S H(t”) t’
dt”}
B(t’) dt’,
(12)
357
where B is the “effective” Finger tensor corresponding to the effective velocity gradient JL?. In the case where ,$ = 0 and H = G, our model thus reduces to Kaye’s [ 111 implicit integral model. Since
trT==Jp2], H and G are functions G, then we have
(13) of tr T. Assuming a linear dependence
2~ T + h,(l+ E’ tr T)T = g,(l
on tr T of H and
+ E’ tr T)Z,
where E’ is a constant. Here to minimize the number of free parameters have chosen the same slope of G and H versus tr T. Equivalently, Xc, ~+(1+E”tr7)7=T-=-4-Q,
2 Go&
(14) we
(15)
where
7
=T-g”Z ho ’
A0 E”
=
E’ 1+3E$’
and
Go
=po (1 -E).
Of course, the free energy cannot be evaluated by eqn. (12) from experimental data of the stress tensor because of the incompressibility conditions. But for the network stress tensor defined in this way, tr r can be taken as a measure of the free energy. Indeed Acierno et al. [ 12,131 have recently proposed a constitutive equation based on this argument. This above constitutive equation is written for a single relaxation time. To take into account for the distribution in N, we allow for multiple relaxation
358
times. The constitutive set of equations #i)
7= E I i
equation
proposed
here is described
by the following
f
Equation (16) represents a spectral decomposition of the stress tensor. Xi and Gi will be shown to be the relaxation spectrum and the rigidity modulus spectrum. They are determined by the equilibrium spectrum, i.e. Q’(W), G’(w) as we will see in later part. The model proposed thus has only 2 adjustable parameters, E and t, both are non-dimensional. 4. Small strain periodic
flow
In a small strain periodic
( )
0
[P]
=+s
flow the effective
1-i
w6 ejwf
velocity
gradient
0.
ejWt
0
0 ‘,
0
0
0,
where 6 is small. Neglect terms of order 0(S2) and write r(p) = r6p) + 6 .jwf ~fp)~ Hence from the constitutive
+
$P
=3
(
1 +
A$ tr
q-LPI
P tb+
Regrouping
+
equation
&
tr
P
S(P)
6
is simply
(16),
&wt
)
7hP)
L,ts2). terms of order 6, it can be easily shown that
(17)
359
The complex viscosity is thus given by (19) The dynamic viscosity and storage modulus are
q’=
c
P=l
G’
=
GPXP 1+
(20)
’
(w&J2
c
wG,h2,
p=l
1+ (CA,)2 *
(21)
It can be readily identified that Xi, Gi are the discrete relaxation and rigidity modulus spectra. 5. Model predictions Figures 2 and 3 summarize the model responses in simple shear and elongational flows for a single relaxation time network. Note that the stress tensor and the time are made non-dimensional with respect to GO and X0, respectively.
The viscometric functions predicted are not sensitive to E values ranging
I
1
2
3
4
I
5
t
6
I
7
I
8
I
9
’ t/b
Fig. 2. Dimensionless viscosity. B = ho7, dimensionless shear rate. - - - E = 0.001, t = 0.2.
e = 0.01, l = 0.2;
360
IT ,~.=lOO
Fig. 3. Dimensionless
G.10
Trouton
viscosity.
G: elongational
rate. E = 0.01,
5 = 0.2.
from 0.001 to 0.01. The onset of non-Newtonian viscosity occurs at about the correct value of the dimensionless shear rate (ho7 = 1). The model predicts stress overshoot qualitatively correctly. Even for high dimensionless shear rates the model does not suffer from the spurious oscillation of the co-rotational models [ 141. In elongational flows, the Trouton viscosity is predicted to be finite for all elongational rates. It is interesting to compare our model response to Acierno et al. [12,13] in transient elongational flows. Ours seems to be better in predicting the stress growth. Theirs predicts a stress overshoot which seems to be structurally impossible in elongational flows. To compare the model predictions with experimental data we have chosen the wide ranging set of data reported [ 51 by the I.U.P.A.C. working party. These data were obtained in six different laboratories participating the program on three types of low density polyethylene, indicated by A, B and C, which had many equal properties including those usually employed for technical characterization. We shall be concerned mainly with sample A at 150°C. It can be shown easily that in steady shear flow the ratio of the second to first normal stress difference is given by N,/N,
= -t/2.
(22)
Due to the lack of second normal stress data for L.D.P.E. we will take t = 0.2 in all subsequent calculations.
361
Fig. 4. Relaxation
spectrum
from dynamic
data taken by Meissner.
The relaxation spectrum, H, is constructed by inverting the function G(t) where G(t) = f 0
H(X) e-‘lh d In X.
(23)
In the following calculations, the relaxation time spectrum will be used in discrete form by taking two time constants Xi in each decade over the range from 10m2 set to lo4 sec. The corresponding rigidity modulus is computed using Gi = H(X,) A In Xi.
(24)
In Fig. 5, the values of the viscosity and the first normal stress difference are reported. The predictions are not sensitive to the values of E ranging from 0.01 to 0.001. The stress growth data and predictions are plotted in Fig. 6 for a shear rate of 1 see-l. We note that the model predicts stress overshoot in shear stress very well. The time at which the first normal stress reaches its maximum is underestimated by the model. Also it is interesting to observe that the model exhibits damped oscillations after the maximum values of the stresses are reached. The results for elongational flows are reported in Fig. 7 at three values of
362
-1
-1
-1
-1
Y IseE’ I Fig. 5. Viscosity and first normal stress difference.
I
0.1
-
MODEL PREDICTION FOR E = 0.001
m
DATA RANGES FROM 8 SAMPLES A-LOPE FROM MEISSNER
1
rs!c, t
Fig. 6. Starting up of a shear flow. 7 = 1 set-l.
100
363
r-
lo9
I
FROM
,
-
DATA
I--
MODEL PREDICTKIN FOR E = 0.01 G = ELONGATIONAL RATE
I
MEISSNER
10' -
/ : :
,---I’
/
I
I
107-
0.1
1
lo
xl0
t (set) Fig. 7. Tensile viscosity.
the elongational rate G in the form of transient Trouton viscosity. The model predicts finite nT for all elongational rates. In fact numerically we find that lim qT = constant X l/c.
G--
(251
Hence the value of E can be optimized using Trouton viscosity data. Such reliable data are not available at the moment for L.D.P.E. and hence we are contented to use E = 0.01. 6. Conclusions We have presented the physical motivations leading to our proposed constitutive equation. This contains two parameters, one is found by the ratio of the second to the first normal stress differences and the other from Trouton viscosity data. The quantitative predictions of the model proposed are good compared to reliable data from the I.U.P.A.C. working party. In conclusion, it seems to us that the network picture of polymeric systems provides a reasonable alternative to the suspension models [ 15-171. The connection between single directive anisotropic fluid and non-affine transformation was shown in a straightforward manner. We have adopted a simple-
364 0.8-
0.6-
B=lO B-12
B=20
Fig. 8. Dimensionless viscosity. E = 0.1, eqn. (26).
minded assumption about the form of g and h. This allowed us to derive the constitutive equations devoid of configuration integrals involving microscopic variables. We have chosen the coefficients of H and G in tr T (eqn. 14) to be the same for simplicity. In fact by choosing different coefficients we can optimize the stress overshoot with regard to the damped oscillation. Indeed by setting the coefficient in G to be zero, i.e. eqn. (12) becomes dT -dt
&T--
.!?TT+
ho(l +E’tr
T)T=g,Z,
(26)
the damped
oscillation is completely smoothed out, see Fig. 8. We believe that one of the main attractions of the model is that there are only two free parameters to be found from experiments. In addition note that one of these parameters (E) may be found from viscometric (“weak”) flows, while the other (e) may be found from an elongational experiment. This seems to be advantageous because then both types of flow response can be controlled.
Acknowledgements We would like to thank Miss A.K. Stewart and Mr K. Frazer for help in preparing the manuscript. RRferences 1 M’. Yamamoto, J. Phys. Sot. Japan, 11 (1956) 413; 12 (1957) 2 A.S. Lodge, Kolloid-Z., 171 (1960) 46; Rheol. .Acta, 7 (1968)
1148; 13 (1958) 379.
1200.
365 3 F.M. Wiegel, Physica, 42 (1969) 156. 4 F.M. Wiegel and F.Th. De Bats, Physica, 43 (1969) 33. 5 J. Meissner, Report of a I.U.P.A.C. working party. To be published in Pure Appl. Chem. 6 D.C. Leigh, Non-linear Continuum Mechanics, McGraw-Hill, Maidenhead, 1968. 7 J.L. Ericksen, Arch. Rat. Mech. Anal., 4 (1959) 231; 9 (1961) 1; Trans. Sot. Rheol., 5 (1961) 23. 8 R.J. Gordon and W.R. Schowalter, Trans. Sot. Rheol., 16 (1972) 79. 9 R.I. Tanner, A.1.Ch.E. J., 22 (1976) 910. 10 AS. Lodge, Mathematics Research Centre, U.S. Army, Technical Report No. 1092. 11 A. Kaye, Brit. J. Appl. Phys., 17 (1966) 803. 12 D. Acierno, F.P. La Mantia, G. Marrucci and G. Titomanlio, J. Non-Newtonian Fluid Mech., 1(1976) 125. 13 D. Acierno, F.P. La Mantia, G. Marrucci, G. Rizzo and G. Titomanlio, J. Non-Newtonian Fluid Mech., 1 (1976) 125. 14 R.B. Bird, 0. Hassager, S.I. Abdel-Khalik, A.1.Ch.E. J., 20 (1974) 1. 15 RI. Tanner and W. Stehrenberger, J. Chem. Phys., 55 (1971) 1958. 16 R.I. Tanner, Trans. Sot. Rheol., 19 (1975) 37; 19 (1975) 557. 17 N. Phan Thien and R.I. Tanner, submitted to Trans. Sot. Rheol.