Constitutive equations for viscoelastic liquids: Formulation, analysis and comparison with data

519

CONSTITUTIVE EQUATIONS FOR VISCOELASTIC LIQUIDS" FORMULATION, ANALYSIS AND COMPARISON WITH DATA A. I. Leonov Department of Polymer Engineering, The University of Akron, Akron, OH 44325 - 0301, USA

1. INTRODUCTION Viscoelastic liquids are capable of accumulating large recoverable strains in flow. This puts the media in an intermediate position between liquids and solids and makes their rheological behavior very complicated. Though the particular properties of elastic liquids have found many useful applications in polymer processing and modern technology, no deep understanding of the nature of viscoelasticity has been reached. Since no fundamental relation is believed to have yet been discovered, at least more than ten popular constitutive equations (CEs) are in competition at present without any clue to a preferable type. For several decades, there have been many attempts to derive CEs for polymer fluids from the viewpoints of mechanics, mathematics and physics. The first approach, which was adopted by many specialists in mechanics, was purely rheological. It postulated nonlinear and quasi-linear relations between the observable variables, the stress tensor cr and the strain rate tensor e. Oldroyd [ 1,2] pioneered the method and also revealed important principles of invariance. The concept was further developed by scientists such as Rivlin, Green and Ericksen and their numerous successors [3]. Later, it was recognized that many rheological equations derived from different approaches are associated with the equations proposed by Oldroyd [ 1,2]. A great many rheological equations, both of differential and integral types have been proposed, and they were able to describe some properties of viscoelastic liquids. The lack of thermodynamic analysis is the main disadvantage of this approach. Hence, such important phenomena as dynamic birefringence, non-isothermal flow, and diffusion

520 cannot be considered simultaneously with rheological constitutive equations. Also, the viscoelastic constitutive equations are often non-evolutionary. The second approach was developed mainly by such scientists as Noll, Coleman, Truesdell, and their colleagues [4-7]. It is the rational mechanics which searched for the most general form of constitutive relations between kinematicand dynamic variables. The basic system involving the constitutive and thermodynamic equations was constructed using strict mathematics, with the CEs having to satisfy the general principles of causality, material objectivity and local action. In this way, the properties of all viscoelastic liquids can be described by a set of hereditary functionals with "fading memory", whose invariance and thermodynamic consistency were perfectly revealed. Unfortunately, there is no unique way to specify the memory functionals and hence predictions are not possible. The third approach is purely physical and explains the behavior of polymeric liquids in terms of the intra- and inter-molecular dynamics. At the beginning, this approach was used to study the behavior of dilute polymer solutions by Kargin, S lonimsky, Kirkwood, Riseman, Rouse, Zimm, Bueche and others (see [8] for a review]). For a long time, concentrated polymer solutions and melts have also been considered as "temporary networks" of entangled chains. Green and Yobolsky [9], Lodge [10], and Yamamoto [1 1] developed semiphenomenological theories that extend the theory for rubber elasticity [12]. This idea was enhanced by the rise of the "reptation" theories, which are due to the work of de Gennes [8], Edwards [13], and Doi and Edwards [14]. In the reptation approach, the creation and decay of the molecular entanglements are studied by the statistical description of a polymer molecule moving along its own axis within a "tube" created by surrounding molecules, and then the motion of the molecule is averaged over high frequency transverse Brownian motion. Another reptational or statistical approach can be found in the works by Curtiss, Bird et al. [ 15-18]. Pokrovsky and Vo|kov [ 19,20] proposed a non-reptational approach based on a generalized Langevin equation for a single macromo|ecule moving in macromolecular environment. A more fundamental approach to the linear viscoelastic properties of polymers has been recently developed by Schweizer [21,22]. This microscopic theory of the dynamics of polymeric liquids called the "mode-mode coupling" (MMC) approach, omits the specific assumptions of reptation theories and formulates the statistical properties of polymers in the most systematic and fundamental way. In contrast to the rational mechanics, the statistical methods have generated a number of specific CEs, and have been generally successful for explaining viscoelastic behavior of polymer liquids in the region of linear or weakly nonlinear deformations. Apart

521

from its poor description of nonlinear data, the physical approaches are also not free from empiricism in formulations and arbitrariness in overcoming mathematical difficulties. Additionally, recent mathematical analyses and numerical simulations revealed examples of unphysical and unstable behavior of the CEs derived by these methods. The fourth method of deriving CEs has been introduced by the author [23,24]. It investigated nonlinear viscoelastic phenomena by the methods of the quasilinear irreversible thermodynamics with the use of a "recoverable strain tensor" as hidden parameter. By this approach, a class of Maxwell-type differential CEs has been proposed under strict stability constraints which are based mainly on thermodynamics. Later, similar approach was adopted by Dashner and Vanarsdale [25] to formulate the general class of CEs almost equivalent to the author's. Almost all CEs proposed in the literature have a limited ability to describe start-up, steady state and relaxation phenomena of polymer fluids in standard rheometric flows, within a relatively narrow region of strain rates usually employed in the tests [26]. However, there are two frustrating problems: (i) none of viscoelastic CE proposed could describe the whole set of available data with one specified set of parameters, (ii) in real modem processing the values of Deborah number, De, may be at least two orders of magnitude higher than those in usual rheological tests, and in that flow region almost all CEs exhibit various instabilities, the reason for which still remains unclear. There are a lot of speculations in the literature on how the instabilities in CEs are related to those observed in flows of polymer fluids There are also contrary opinions about the physical sense of non-evolutionary behavior of viscoelastic constitutive equations (see, e.g. Refs. [27,28] with contrary views on the subject). The temptation to relate the instability in CEs to real flow instability was perhaps caused by the widely spread perception that no rheological constitutive model is globally stable. A few books [ 18,26,29-31 ] and a lot of papers on the polymer rheology do not answer the question as to which CE should be chosen to solve fluid mechanical problems of polymer processing in the usual case of large De number flows, where the nonlinear effects of elasticity are important. Therefore, general studies in this field should be aimed on searching for not only descriptive but also reliable CEs. Thus in our opinion, the principles of choice of CEs should also include, besides the descriptive ability of viscoelastic CEs, such fundamental thermodynamic properties of these as formulation of dissipation and free energy, along with stability constraints.

522 To study these fundamental properties, even for the relatively simple CEs, one has to employ a general formalism within the framework of which it is feasible to find some general stability constraints. Several formal approaches to the viscoelastic CEs have been proposed relatively recently in the literature. One of them is the local quasi-linear approach of non-equilibrium thermodynamics developed by the author [23,24,32]. As mentioned, this approach has resulted in derivation of a class of Maxwell-like CEs, with some thermodynamic constraints. Almost 10 years later, the Poisson-bracket variational approach has also been established. It was first introduced by Grmela [33,34] and then elaborated by Beris and Edwards [35]. This approach, extending the Hamiltonian formalism in classical mechanics to the case of continuum mechanics, also employs the dissipative functional. It was proved [36] that both the above approaches result in the same "canonical" formulation of general Maxwell-like CEs. Another two general formal schemes, consistent with the irreversible thermodynamics, have been developed by Jongschaap et a1.[37,38] and by Kwon and Shen [39]. The first one employs a matrix representation in the phenomenological relations between the thermodynamic forces and fluxes, with additional relations between external and internal (hidden) thermodynamic variables. The second one uses the notion of evolution of the temporary network structure. Regardless of all the differences in their detailed schemes, their equivalence to the general author's approach is evident. Another objective that should be kept in mind when proposing the CEs for polymeric liquids, is solving geometrically complicated flow problems under high deformation rate, as related to industrial processing. Complexity of flow problems for polymeric liquids is enhanced by the effects of fading memory which do not exist in viscous liquids or elastic solids. Hence, viscoelastic polymer fluids demonstrate many unique rheological effects, such as nonlinear evolution of stress under steady deformation, which cannot be seen in the other fluids. Because of this, even geometrically simple problems which may be solved analytically for viscous fluids, have often to be treated numerically in the viscoelastic case. This necessitates the elaboration of certain criteria for the selection of CEs for practical applications. The following principles of choice have been suggested [40,41 ] for the selection of CEs: (i) Stability. However well an unstable CE can describe rheometric tests, it is impossible to use it in modeling of polymer processing, since the Deborah numbers there may be at least two order of magnitude higher and the flow much more complicated. Extrapolation of the majority of CEs to the region of high Deborah numbers and 3D flows usually results in several types of instabilities in numerical flow simulations. These instabilities reflect the mathematical

523

structure of the proposed CEs. In the most cases, they are not related to the physical instabilities observed in the flows of polymeric fluids, or poor numerical algorithms, but rather to violations of some fundamental principles. (ii) Descriptive ability and flexibility. It is now well recognized that polymer melts with similar linear viscoelastic spectra can show qualitatively different nonlinear behavior. For the proper description of various flows, this requires some functions of the kinematic variables, and the associated nonlinear parameters in the CE (which vanish in the linear limit) to be specified within the stability constraints. Once these functional forms and parameters are specified for a particular polymer, the CE must simultaneously describe the entire set of available experimental data fairly accurately. (iii) Computational economy. The proposed CE should allow for numerical calculations in complex flows with as little computational effort as possible. For example, despite the good descriptive ability, it is rather cumbersome to work with models in which the elastic potential is specified in terms of the principal values of a strain measure, and it is usually conceived that working with CEs of differential type is preferable for numerical calculation than with integral ones. (iv) Extensibility. Real polymer processing is confronted with a variety of complications such as compressibility, non-isothermality, wall slip, phase transitions and separations, chemical effects (degradation, curing), etc. In principle, the CE of choice should be amenable to extension in order to accommodate these phenomena. Among the four principles listed above, the first two can be regarded as the most fundamental properties which the CEs should possess.

2.

FORMULATION EQUATIONS

OF

VISCOELASTIC

CONSTITUTIVE

There are two types of viscoelastic constitutive equations which are only in use today for practical applications. These are of differential (Maxwell-like) and single integral types. Therefore we discuss below only the formulation of these two types of CEs, since only they have been tested up to present.

2.1. Formulation of general Maxwell like constitutive equations [36] 2.1.1. Thermodynamics and general single-mode viscoelastic approach We will employ in this section the general approach of irreversible thermodynamics discussed in various books (see, e.g. [42,43]). The fundamental hypothesis of"local equilibrium" underlying the method assumes that even in a

524 non-equilibrium process, in any arbitrarily small macroscopic particle of a medium, it is feasible to operate with the same common thermodynamic functions, depending on the same variables, as in true thermodynamic equilibrium. It is also possible to involve in the analysis some non-equilibrium parameters which vanish in equilibrium; the less these are involved the more economic is the description. This assumes that the system under study is in a sense close to a thermodynamic equilibrium. For viscoelastic liquids, the true thermodynamic equilibrium is at the rest state. The same holds for elastic solids, but for the liquids there additiona|ly exists an incomplete thermodynamic equilibrium in the stressed state. Though strong flows of elastic liquids can be far away from true thermodynamic equilibrium, we can still assume that they are close to the state of incomplete equilibrium which is characteristic for stressed elastic solids. In fact, this fundamental hypothesis underlies all the theoretical treatments of elastic liquids, the molecular approaches, Poissonbracket and matrix formalisms included, and reflects our intuitive view of shear elasticity in liquids. We assume that in the simple case under study, the state variables of elastic liquids are" the temperature T and a hidden variable, symmetric second-rank non-dimensional tensor =c. We assume additionally that the tensor ___eis positive definite, which to some extent, have been justified independently [36,44]. In the following, we employ the Eulerian formulation of the constitutive equations and, without loss of generality, use a Cartesian coordinate system. The physical sense of tensor c can vary, however. In our own studies [23,24] the tensor =ewas treated as the Finger elastic (recoverable) strain C, 1, which may be measured independently. When related to molecular approaches, the tensor c is proportional to the averaged diadic built up by the end-to-end vector of a part of a macromolecule between two entanglements. In these approaches, the configuration tensor is a typical internal variable and the problem of how to measure it directly is questionable. Being mostly interested in studies of isothermal flows of isotropic elastic liquids, we can now introduce as a proper thermodynamic potential, the free energy density per mass unit, F = F( T, It, 12, 13). Here /1 - tre,

12 -

1/2(I12-tr__e2),

13 - detc,

(1)

are the basic invariants of the tensor g. Then the "thermodynamic stress tensor" 0% associated with the free energy F, as in the quasi-equilibrium elastic case [45], is introduced as follows" (o~ = 2pg=.c3F/Oc = 2p[F,c_- + F2(I,c=- s + Ffl,6=] (Fj = OF/OIj) (2)

525

Here 9 is the density and fi is the unit tensor. The only physical reason to operate with the equation (2) is the explicit assumption [23,24] that polymeric liquids always possess an "elastic limit", a quasi-equilibrium situation achieved on very rapid (instantaneous) deformations, where the temporary entanglements in macromolecules act like cross-links in cross-linked rubbers. Thus, in this limit, when g --+ _Cl the use of equation (2) can be justified. It should also be noted that no assumption of incompressibility has been made. We now briefly discuss the dissipative effects [23,24] which are associated with flows of elastic liquids. By using the laws of mass and energy conservation as well as the momentum balance, and separating the total entropy variation in the flux and origin ("entropy production"), one can easily obtain the ClausiusDuhem expression for the entropy production Ps in the case under study: TPs = -~I.VT + tr(=o-.e) - pdF/dt [ T

(3)

Here ~/is the thermal flux, __o-isthe actual stress tensor, g is the strain-rate tensor, and the symbol d/dt means the time derivative in the frame of reference associated with a moving particle of liquid. The first term on the right-hand side of equation (3) reflects non-isothermal effects and could be considered as independent of the other two, which are related to mechanical dissipation. According to the Second Law, P~ should be positive for all non-equilibrium processes and vanish at equilibrium. Considering further the isothermal behavior of elastic liquids. Using equation (2), we can rewrite the dissipative terms in (3) as follows" D = TPs[T = tr(cr-e) - tr(_O-e'e-l"1/2d__e/dt),

(4)

where D is the mechanical dissipation. The next step in developing the non-equilibrium thermodynamic approach is to represent the dissipation as a characteristic bilinear form D = ~ k " Yk where Xk are thermodynamic "forces" and Yk are conjugated thermodynamic "fluxes". As thermodynamic forces, we can naturally take the actual, o-, and equilibrium, =O-e, stresses. Then the thermodynamic flux conjugated to the actual stress =o-is certainly the strain-rate tensor g. Another thermodynamic flux, ee, conjugated to the quantity ere is a tensor which is related but is by no means equal to 1

-c 2 =

-I

.de/dt.

This is because the latter does not even satisfy the frame-

=

invariance conditions. This situation demonstrates an ambiguity in definition of

526 the second thermodynamic flux, e~. Mathematically, it associates with arbitrariness in finding a tensor from a scalar product of two tensors. Thus an additional physical assumption is needed to define the unknown quantity e~. This has been proposed in [23] as follows: not only thermodynamic force, cr~ but also conjugated to this, the thermodynamic flux ee, have to be defined as in the case o f quasi-equilibrium elastic solids. This assumption seems to match reasonably well the local equilibrium hypothesis. It also results in the following o

procedure" to find the solution, e~, from the kinematic relation c=c.__e+e.__c,valid for elastic solids (or for the total continuum) and hold the expression for the o

non-equilibrium case. Here c is the co-rotational or Jaumann time derivative of the tensor c. It is exactly the same procedure that we have used in the definition of the equilibrium stress tensor ere. In so doing we can define the quantity e, as the solution of the equation" o

f = C.ee + ee'_C_,

(5)

which has been obtained in an explicit form [24], too awkward to be reproduced here. It is easy to see from equation (5) that the term pdF/dt I T in equation (3) is reduced to tr(o-.e~), which, in turn, can be represented in the form of the second term in the dissipation equality (4). It should also be noted that this derivation, being general enough, was confirmed by an independent specific analysis of kinematics of viscoelastic deformations [24,46,47] when ___e= C~l. If the nonequilibrium stress _O-pand strain rate ep are defined as follows, _O-p- g - o-~,

__ep- __e-ee

the dissipation inequality due to equation (5) is rewritten in the form: D - TPs[T - tr(_O-p.e) + tr(o-~.ep).

(6)

Equation (6) is now represented as the typical bilinear form D = f~(k'Yk discussed above, where C~pand ee are independent thermodynamic forces, and ___ep and ee are independent thermodynamic fluxes. Two independent sources of mechanical dissipation are now clearly seen from equation (6); both of them being positive defined and vanishing in the equilibrium. They are" (i) the power produced by the irreversible stress _O-pon the total strain rate =e (the first term in (6)), and (ii) the power produced by the reversible stress o-~ on the irreversible

527

strain rate ep (the second term in (6)). Using quasi-linear scheme of irreversible thermodynamics, these were connected in [23,24,46] by phenomenological relations with kinetic coefficients represented as some rank four tensors depending on the internal parameter, tensor g. In presenting these relations, the Onsager symmetry of kinetic coefficients, proved for quasi-linear case in [43], has also been used. The structure of the kinetic rank-four tensor had also been completely revealed [24]. An interested reader can find the details in papers [24,46,47]. This strict, straightforward and complete approach can be compared with somewhat controversial matrix approach [43,44].

2.1.2. Maxwell models with quasi-equilibrium stress. We will consider from now on only the Maxwell liquids for which the actual stress tensor is equal to the thermodynamic stress, i.e. =o-= fie- (=O-p= 0). In this case, the first term in equation (6) vanishes. Then the quasi-linear phenomenological relation between the thermodynamic force, _fie and the generalized thermodynamic flux, gp, is of the form" ep,ij = m0kl(T,c)O'e,kl(Y,c).

(7)

Here M is a rank-four mobility tensor, an isotropic tensor function depending on the configuration tensor =e; the stress tensor _o-~now having the form shown in equation (2). Substituting (7) into the dissipative equality (6) with O-p = 0 represents the latter in the quadratic form of the stress tensor, D = M~st-cr~-=o-~.~t. This in turn results in the fact that the mobility tensor has to be positive definite and has the following symmetry properties" it is symmetric in the first two and second two indices and in transposition of these indices. Also, substituting expression (7) for ge into equation (4) yields the evolution equation for the configuration tensor c" V

c+~c)/0(T) = O,

~c)/O(T) = 2C-ep(T,c).

(~j(T,c)/O(T)

=

2CikMkist(T,c)O's,(_~._C))

(8) V

Here c is the upper convected time derivative and ~ c ) is a non-dimensional isotropic tensor function of tensor c defined as in equation (8) and related to the dissipative processes in an elastic liquid. Additionally, in the rest state when c--+ =_8(I1--+I2---~3; I3---~1), we have" =d(T,c)---~=0. We assume that this limit transition is regular, meaning that there is a limit to the linear Maxwell viscoelastic equation when the intensity of strain rate is very low. In the incompressible case with

528

linear dependencies Cre(=e)and =d(=c), equations (2) and (8) are easily reduced to the upper convected Maxwell model. Multiplying the first eqation (8) from the right and from the left by g-~ reduces the latter to the a "dual form'" + s

= 0,

(b = __e-')

(8a)

where b is the lower convected time derivative of the tensor b. Also, in the incompressible ease with linear dependencies of ~ and ~ ) , eqations (2) and (8) are reduced to the usual form of the lower convected Maxwell model. The difference between these is mainly in the dependence of stress on elastic strain" =o-= G=e corresponds to the case of rubber elasticity, whereas =o-= -Gb can describe the elasticity of crystals or metals. Equations (1),(2) and (8) give rise to the energy relation:

r

+ D - tr(Cre.e),

D - tr(g~- Cre.gb)/(20(T)),

(9)

where D is the mechanical dissipation. If the new kinetic rank-four tensor L(e=) is introduced as follows" Lijkl = CimCjnMnmkl

(10)

the dissipative term ~___e)in equation (8) is rewritten in the form"

~j( e=)/O(T) = LoglpOF/c3ckl

(11)

Substituting equation (11) into (8) represents the latter in the form which coincides with the evolution equation obtained by Beris and Edwards (see equation (2.13) in the second part of Ref.[35]). The tensor L has the structure [35] similar to that known for the tensor M. We now analyze the compressibility condition in the general Maxwell model described by equations (2) and (8). Multiplying equation (8) scalarly by (2=e)1 gives:

1/2d/dt(lnl3) + trep - tr__e; trep - 1/2tr[__c-l.=d(_c_)/0(T)], (12) where the function =d(c) was defined in equation (8). There is also the mass conservation equation which can be written in the form"

529

d/dt(lnp/po) = -tre.

(13)

Here/90 is the density in the rest state. Then, combining equations (12) and (13) yields:

d/dt[ln(,O/po. ~-3 ) + trep = O.

(14)

Equation (14) shows that P---~Po when c--+_6, since in this limit, trep--+0. It also shows [23,24] that if

trep = 0,

#/9o = 1/-4~-3,

(15)

exactly as in the equilibrium limit of elastic solids. Thus, the evolution equation (8) for tensor ___calso describes the law of mass conservation in compressible elastic liquid only if the condition tr__ep = 0 is satisfied. When trep ~ 0, the density variations are not described anymore by the configuration tensor c but satisfy the kinetic equation (13) or the equation of mass conservation (12). It means that, in this case, the density is not a state (thermodynamic) variable and an attempt to improve the situation by making the assumption that F = F(T,p,e=) contradicts the Murnaghan's definition (2) of the thermodynamic stress. This was the only reason why in our publications [23,24,46] we also employed the condition trep = 0. A simple constitutive equation for the compressible part of stress has been discussed in [23,24]. It includes the equilibrium pressure and also KelvinVoight modeling of compressible effects. This approach has been recently extended for non-isothermal linear compressibility phenomena with a bulk relaxation spectrum when introducing a set of scalar independent hidden parameters {~k} [48]. Consider now the incompressible case when p =/9o - const. When tr__ep :/: 0, equations (8) still hold, whereas equations (13) and (14) are reduced to: 1/2d/dt(ln/3) + tr__ep- O; tr__e= O.

(16)

Surprisingly enough, in this case, the free energy remains the same, i.e.

F=F(T,II,I2,I3), because the tensor =e is not directly related to the density. By introducing the modified free energy / ? - F - p represented by

lnI 3, the actual stress =0-can be

530 __o-=-p=~+ ere.

(17)

Here the isotropic term 'p'(a Lagrange multiplier) in equation (17) is introduced to satisfy the condition of incompressibility tre = 0; the thermodynamic tensor ere being still defined by equation (2). This means that the tensor cr~ is not defined with the accuracy of an isotropic term as in the equilibrium case, and isotropic pressure is not an equilibrium one. These features result in the fact that, even in the incompressible case, the expression (6) for the dissipation is now not invariant under the transformation, _o-e--~o-~+p__6. When the constraint trep = 0 is used, equation (15) results in the incompressibility condition [23,24],/3 = detc = 1. In this case, the modified elastic potential is introduced as" W = W(T, II,I2)-f~(I 3- 1), where W = pie. The ^

magnitude ~ here is the isotropic pressure. It is used in differentiation in equation (2) as the Lagrange multiplier, allowing us to consider all cij variables as independent. Then equation (2) results in the Finger formula for the stress: - -p6_+ 2 Wle- 2 W2c -l.

(18)

Here Wk= OWlOIk and p = ~-/;'2/2. Also in this case,/2 = trc l and the formulae for free energy and the stress tensor are exactly the same as for elastic solids. Also, isotropic pressure and the expression for dissipation (6) are invariant under the transformation: _O-e--)Cr~+ p=6. We finally stress that as shown in this section, any attempt to extend the configuration tensor approach to the compressible case is thermodynamically inconsistent. Since the Maxwell type models involving configuration tensor c were derived from the evolution equation for distribution function, it simply means that the later equation should be somehow modified. It is unknown at present what should be done to rectify the situation.

2.1.3. Maxwell models with non-equilibrium and non-potential stress tensor[36]. The above thermodynamic derivation of Maxwell-like constitutive equations was further extended [36] to include formally into consideration the incompressible Gordon/Schowalter [49], Johnson/Segalman [50] and PhanThien/Tanner [51,52] models. For this reason, instead of the natural, quasiequilibrium relation (5), the following non-equilibrium evolution equation for the tensor c was proposed"

531 o

c = ~(__C'ee+ ee'__c)

(5a)

Here ~ is a numerical parameter (-1 <~< 1, %~0). The case %= 0, related to the corotational derivative, can be also easily considered but we shall omit this. Then following the same procedure, described in the derivation of equation (8), we can finally obtain [36]: o

__o -= -p__6+ __ere• Crex= Cre(g)/~,

_c- ~(e.e + =e.___c)+ &(T,g)/O(T) = 0,

(19)

D = tr(Crex-&). In equations (19), written here in so-called canonical form, the formula for the extra stress tensor Crex is given due to some thermodynamic reasons [26]. Note that when ~,~+1, equations (19) have no perfect elastic limit and therefore they are non-equilibrium. One can also ignore the thermodynamic origin in the stress formulation ere(C) and consider formally the second formula in equation (19) as a general nonpotential stress-elastic strain relation O-~x= Crex(C,~). It can be formally done even in the case ~,=+1, where there is a limit to elastic solid behavior, ~(T,_e_)--+=0. However, in this case, it is always possible to consider such a loading-unloading procedure that will create work from nothing [53]. It means that this approach inevitably leads to the theoretical eternal motion machine and therefore is thermodynamically forbidden.

2.1.4. Examples of single mode Maxwell-type constitutive equations. We now illustrate how the particular Maxwell-type CEs can be obtained from the above general equations by specifying the terms ~ and ere in equations (19) and (2), or F in equation (1). When ~b~= g-_6,

(2p/G)F = I~-3,

ere = G__e,

(20)

where p and G are the shear modulus and density, equations (19) represent an interpolated Maxwell or the Gordon/Schowalter [49] - Johnson/Segalman [50] model. When ~=+1 or ~=0, it is respectively called the upper/lower convected or co-rotational Maxwell model. When

532

0 = 0o(T)f(I1),

g~ = g-=b,

(2p/G)F = I~-3,

_O-e=G___c,

(21 )

it becomes the general Phan-Thien/Tanner model [51,52], where f(I1) was proposed as a linear or exponential, increasing function of I1,. Again, when ~=1, we particularly call it an upper convected Phan-Thien/Tanner model. When 0 = Oo(T)f(IIe), IIe=tr(__e2), ~= 1, ~=__c-__6",(2,o/G)F=Ii-3, O-e=Gg,

(22)

where f is a decreasing function of lie, it reduces to the White-Metzner CE [54], which does not belong to the general class of quasi-linear CEs. The FENE model (see, e.g., Section 8.5.3 in Ref. [26]) can be written as: Y_,=1, 0==K(I~)c-=d, (2p/G)F=(Ic3) 1nK(I,), o-~=GK(I~)c, K(I, ) = (Ic- 3 )/([c-ll ),

(23)

where Ic=Rcz and Re is the ultimate finite dumbbell length. When ~=1, =dr=otc2+(1-2ot)c-(1-c~)__d(0
__o-~=Gs

(24)

it represents the Giesekus model [55], with a numerical parameter or. The case ~=1, =d~=B(I~)(c-__6), (2p/G)F=(3/~)lnB(I1), B(I1)--1+~(I,-3)/3

_O-e=G__c/B(I1), (25)

corresponds to the canonical representation [36] of the Larson differential model [56]. The general class of the author's incompressible CEs [23,24,47] is represented as:

~=1, r l /2c[b~(c-I~ 5/3 )-b2(c~-12~_/3)], O-e=2Wlc-2W2c-~

(26)

Here W(I~,Iz) (=pF) is the strain energy function for incompressible case. The positive functions bi(Ii,I2) should have a proper linear viscoelastic limit and their positive definiteness suffices the positive definiteness of the dissipation. The convexity constraints,

533 Fl>0,

F11Fzz>F12 z, (Fi=c3FIc3Ii,Fij=c3Filc3Ij)

F2>0,

(27)

imposed on the general form of potential F were also suggested [23,24,47]. The important implication of inequalities (27) and the proper use of this class of CEs are discussed in detail in [40]. In the simple case of bl=b2=l, with the neoHookean potential for F, it reduces to the simplest "Leonov model" which includes no nonlinear parameters.

2.1.5. Multi-mode Maxwell models. All the above single mode viscoelastic constitutive equations of Maxwell-type are usually extended to the multi-mode case. We briefly discuss below only models with potential stresses in the common incompressible case. The typical N-mode extension is as follows" N

F - y'Fk(T, ck), k=l

N

cy - ~ r ---ex

-=

(T,c ) --- ex,k

=k

(28)

~

Here Fk, Crex,kand Ck are the free energy, extra stress tensor and configuration tensor, respectively, in "k"th relaxation mode. Additionally, the general form of evolution equation presented in equation (18) holds for every configuration tensor s In the limit of very low Deborah number, this approach shows the linear viscoelastic behavior, with a discrete spectrum of relaxation times { Ok}. It means that every nonlinear relaxation mode is generated by the corresponding linear viscoelastic mode. The generalization to the multi-mode case can be justified only if the various relaxation modes are well separated, i.e. 01>>02>>...>>ON

(29)

In this case, it is reasonable to assume that the various relaxation modes act independently. All the known experimental datatestify in favor of inequalities (29). Additionally, the guess-independent Pade-Laplace method (see e.g. [57]) reveals the effective discrete linear relaxation spectrum in accord with (29).

2.2. Formulation of nonlinear single integral constitutive equations From a wide class of viscoelastic CEs of the integral type, only the single integral ones have been experimentally tested. In the common incompressible case, its general form is represented as [58]"

534 !

o%•= ~ [qg,(I,,I2,t- x)C- qg2(II,I2,t-

"r.)C-l]dx

(30)

--00

Here _O-ex is the extra stress tensor, C is the Finger total strain tensor for incompressible media, whose time evolution is described as follows" v

C-dC/dt-C.Vv-(Vv_) _ _

_._

_~_

r -C-O; ~

___

=

C[

=8 I =

T

(31)

- -

Here I l= trC, I2 = trC -1, and ~, and q32 are generally independent functions. We can also introduce many other measures of deformations. One of them, the Hencky measure, H = (l/2)lnC, will be used below. Experiments show that the simplified time-strain separable version of equation (31),

6k (Ii,/2,t-x) = m(t-x)q)k(I1,/2)

(k=l,2)

(32)

can be introduced. Here re(t-x)- dG(t-x)/dx, and G(t) is the relaxation modulus. Equations (30) and (31) are not the only single integral form of CEs. When, for example, the mixed convected time derivative is used, the CE can also be represented in an integral form (see, e.g., Ref [26]). Also, Kaye [57] and Bernstein et al. [60] proposed the potential form of equations (30) and (31):

(O, =(2p/G)OP/DI~,

~p2=(2,o/G)0~'/012;

(33)

or in the time-strain separable case"

q)~=(2p/a)3F/Ol~,

tp2=(2p/a)OF/OI2.

(34)

Here G is the elastic Hookean modulus. Equations (33) or (34) constitute the KBKZ class of single integral CEs. The potential F in equations (33) denotes the thermodynamic free energy with relaxation effects taken into account. For the time-strain separable viscoelastic CEs with potential F, the basic functionals such as the free energy W, the extra stress tensor o-e and the dissipation D are of the form [61 ]: ^

535

(35)

W = (p/G) iF ( I, ,I 2 )m(t - z)dx -oo

t

(36)

o-~ = (2p/G) IC. c3F/~C(I,,I 2 ) m ( t - , ) d , --O0

t

D - tr(crce)-dW/dt = Co~G) IF (I~,I 2 ) [ d m ( t - x) / dxld~

(37)

--oo

A comparison of equations (2) and (8) with equations (35) and (36) clearly shows that CEs of the differential type where the dissipation and free energy are generally independent, are more flexible for rheological modeling than CEs of the integral type where those quantities are roughly proportional. We now describe some particular CEs of separable single integral type by specifying q~, and q~2 in equation (30) or the potential F for the K-BKZ type in equations (34). Wagner et al. [62] proposed their first specification as: q~l=fexp(-nl 4I - 3 )+(1-J)exp(-n2 ~/I - 3 ),

I=flll+(1- /2;

r (38)

where f, nl, /
q 2:

h(I~,/2) = l/(ot-,/~), z = (/1-3)(I2-3). (39)

Here again, a and, , are positive fitting parameters. We further refer to the specification (38) as the Wagner model I, and (39) as the Wagner model II. Luo and Tanner [64] proposed the modification of the CE originally presented by Papanastasiou et al. [65] in the form: qo,= oth(Ii,I2),

qo2= ~:oth(IiJ2);

h(IiJ2) = [ot-3+flll+(1-tS')I2]1,

(40)

where c~, /3 and ~: are fitting parameters and K>0. When K=0, this model becomes identical to the one proposed by Papanastasiou et al. [65]. One can easily see that no thermodynamic potential relation exists (q012 ~ q~21) in all of the above three CEs (38)-(40).

536

The Lodge model [10], which is reduced to the integral presentation of UCM[1], q01=l,

q)2=0,

(41)

is a particular case of the separable K-BKZ class (34) when the neo-Hookean potential is used. In order to describe better the viscometric data, Larson and Monroe [66] suggested the following 4-parametric form ofF for the K-BKZ class: (2p/G)F = (3/2ct)1 n[ 1+~(I-3)/3], I = (1-fl)ll + ~ / 1 + 21312 - 1, Ot=Ko+K2tanq[K1A3/(1+A2)], A=I2-I1

(42)

Here Ko, K1, ~:2 and fl are numerical fitting constants. Another potential form was derived by Currie [67] as a close approximation of the Doi-Edwards reptation model [ 14]. It is written as (2p/G)F = (5/2)1 n[(J- 1)/7],

J = ll + 2x/I2+ 1 3 / 4 ,

(43)

and contains no nonlinear fitting parameters. A linear combination of simple potential forms, t G2 e x p ( - t 2p/~=. G! exp( -~- ). In[1 + ~( I , - 3)] + 21302 -~2 )" In[1 + 13(12 - 3)] + 2or0 -

G3 exp( - t 20--7 O3 )" ( I , - 3)

(44)

was also introduced by Yen and McIntire [68] as a partially time-strain separable version of the general potential presented by Zapas [69].Here the relation between the potential /~ and CE is shown in equations (33), Gi's and 0,'s are moduli and relaxation times, respectively, and a and fl are nonlinear fitting parameters. It should also be noted that the thermodynamic potentials proposed in the theory of rubber elasticity can in principle be applied in the case of K-BKZ class of integral CEs as well as for CEs of the differential type. In many cases, they are presented as a function of principal values Ci of the Finger tensor rather

537

than invariants of the total Finger deformation tensor. Below are several examples of these potentials.

(i) Ogden potential [70]: 2PF= ~ G,, (C.,,~ 1 2 +C~",~2 + C ~,,~2 3 _ 3)

(45)

Here Gtn are numerical parameters which can be negative or positive, and Gn are shear moduli. The potential (45) becomes identical to the Mooney potential [71] if it contains only two modes corresponding to a 1=2 and cx2= -2.

(ii) Valanis/Landel potential [72]: 3

2pF= i~[aff-Q, (lnff-C , - 1 ) + 131n.fQ-~]

(46)

where c~ and 13 are parameters with dimensionality of modulus.

(iii) The BSTpotential [73 ] 2pF-- A I ~ + B _

_

n

o

m IE,

I E - 1 (C~,,/2 , + C ,~2 2 + C ,3/ 2 ) . ~

n

n

(47)

n

Here A and B are the parameters with the dimensionality of modulus, whereas n and m are numerical parameters. As mentioned, the non-potential viscoelastic CEs proposed in the literature, are unphysical, since when applying very fast deformations, it is possible to create a perpetual motion machine from a hypothetical material subordinated to this type of CE [53]. Nevertheless, we will analyze below this type of CEs too, considering it as mathematical abstraction for the superficial data curve fitting. To make the following analyses more efficient, we now introduce a unified set of notations for both differential and single integral types of CEs which employ only upper convected time derivatives in the evolution equations. The lower convected time derivative (~=-1) in equation (19), can be equivalently rewritten in the form of upper one (~=1) [36], as shown by equation (8a). We introduce a modified pressure term defined as:

538

for differential CEs P

t

+ jq)212m(t- 1:)d1:

for integral CEs

(48)

--cO

Then using the Cayley-Hamilton identity and the invariance of rheological variables (say, the extra stress o-~) under arbitrary addition of isotropic terms in the case of incompressibility, we represents both the classes of CEs as follows: ____o-=-p'=d+ O-e,

O-~= 2,~__.SF/SR, for differential CEs

o

=e

[_~m(t_,)E(c)dx

for integral CEs

(49)

E = 2p[q01c + q02(Ilc - c 2) + q0113~ Here 8/8g is in general the partial Frechet derivative with respect to c, and with the definition for q)i in equations (2) and (30). In this notation, c becomes the total Finger strain tensor __Cin the case of integral CEs, and thus for the CEs of integral type and the author's CEs, I3=1, q)3=0 automatically. Even though the set (40) is written for hyper-viscoelastic equations, the non-potential viscoelastic formulations can also be included into consideration. The (elastic or total) strain tensor c is the solution of the following evolution problem: t el,__o = =8

v

for differential CEs (50)

= = =

=

[ Ct=t,

= C] t = t t = =8,

~(c) - =0 __ =

for integral CEs

Here ~bis the dissipative term which vanishes for integral CEs, and from now on we consider only upper convected time derivatives even for differential models.

3. STABILITY OF VISCOELASTIC CONSTITUTIVE EQUATIONS

3.1. Background In numerical simulations of viscoelastic flows, degradation of the numerical solution or lack of convergence of computational schemes has been frequently

539

observed for large or even modest values of Deborah numbers. It is thought that the main reason for this numerical malady is an improper choice of a CE which possess some bad mathematical properties (e.g., see Ref. [74], p.314). Therefore the analysis of such properties of various viscoelastic CEs as their stability and boundedness seems to be a pre-requisite for any successful numerical modeling. It will be shown in this Section that these properties are in turn, originated in our capability to incorporate basic thermodynamics laws in formulation of the CEs. We discuss in this section recent results of stability analyses for isothermal formulations of CEs for viscoelastic liquids. To begin with, we briefly describe the basic aspects of stability analysis in general, and two types of instabilities related to the formulation of viscoelastic CEs in particular. The general purpose of the stability analyses is to find the conditions of boundedness for a solution of a set of functional equations subordinate to some additional (usually initial and boundary) conditions. In our case, the total set of equations consists of momentum balance and continuity equations coupled with viscoelastic CEs. No-slip boundary conditions are usually applied at the rigid surfaces and the common set of kinematic and dynamic boundary conditions is used at the unknown free surfaces. It is assumed a "basic" solution of the set (generally, 3D and time dependent) to exist, subordinate to initial and boundary conditions. The solution is said to be stable in the Lyapunov sense, if "small" disturbances imposed on the solution at time to, which satisfy the homogeneous boundary conditions, remain to be small at any time t (to
540 The purpose of the stability analyses discussed in the following sections is not to describe the real physical instabilities observed in flows of viscoelastic liquids, but rather to reveal some constraints which should be imposed on the formulation of CEs to prevent them from the occurrence of unphysical instabilities. Therefore, we call the stability considered in the subsequent analyses not the stability of viscoelastic flows but the stability of CEs. These are related to the formulation of CEs. The formulation of nonlinear viscoelastic CEs is not as straightforward as it appears at first sight. For many specifications of CEs in the above two classes, numerical modeling of high Deborah number flows has displayed heavy unphysical instabilities, regardless of the "true physics" employed in their derivation. Two types of instabilities related to the formulation of CEs were observed for large or even modest values of Deborah numbers, and have been analyzed in the literature. These are the Hadamard and dissipative instabilities. The case of Hadamard stability is the most understood. We define the complete set of equations for viscoelastic liquids as Hadamard stable (or evolutionary, or well-posed) when the solution of the boundary-value Cauchy problems for the set at any time provides the complete initial conditions for determining the solution at subsequent instants in time [75]. Thus, the Hadamard stability allows one to continue the solutions in the positive direction of the time axis. When this is impossible, very quick blow-up instability occurs, with extremely short wave disturbances, which results in progressive failure in numerical calculations: the finer the mesh, the worse will be the degradation of the results [76]. In many cases, one can treat the Hadamard instability as a blow-up type increase in the amplitude of initially infinitesimal waves of disturbances as the wavelength tends to zero. In viscoelastic liquids, this type of instability can be associated with a nonlinear rapid response of CEs. Hence, this type of instability depends on such quasi-equilibrium properties of the CEs as the type of differential operator in the evolution equation for differential models and the elastic potential in the hyper-viscoelastic case. Rutkevich [77] initiated the study of Hadamard stability for viscoelastic CEs. Later, Godunov [78] analyzed some aspects of the stability in more detail. Some significant results were obtained relatively recently by Dupret and Marchal [75] and Joseph and co-workers (see, e.g., Joseph's monograph [76]). These studies which have analyzed the Hadamard instabilities in particular flows for particular specifications of differential viscoelastic CEs, can be summarized as follows. Joseph and co-workers [79,80] and independently Dupret and Marchal [75] proved that the interpolated Maxwell model, which involves the mixed time derivative in the evolution equation, is Hadamard unstable, except for the

541 cases of upper and lower convected time derivatives. The ill-posedness of the Gordon/Schowalter [49] or Johnson/Segalman [50], and the original PhanThien/Tanner [51,52] CEs is subject to this type of instability. Using the general method of characteristics, Dupret and Marchal [75] also showed that the WhiteMetzner model [54] is nonevolutionary, which was again justified by Verdier and Joseph [79], who employed a perturbation method in their analysis. It was found that the dependence of the relaxation time on the invariant of strain rate tensor is the cause of Hadamard instability in the White-Metzner model. Analyzing the White-Metzner model in elongational flow, Verdier and Joseph [81] also noticed a type of dissipative instability that occurs whenever the extensional strain rate exceeds the half of the reciprocal relaxation time. The method of characteristics is the most general mathematical tool for the Hadamard stability analysis. This method is, however, sometimes complicated and cumbersome. Fortunately, in the most interesting cases, it can be simplified by using the "frozen coefficient" method (see, e.g. [76]), when analyzing only extremely short and high frequency wave disturbances propagating with a finite speed. Those cases are related to all the CEs of the quasi-linear differential type as well as the time-strain separable single integral CEs. Then, the linear stability analysis of the problem can be studied locally, without considering boundary conditions. Although following Kreiss' examples [82], this local stability condition is neither necessary nor sufficient for the overall stability in the general case of nonlinear partial differential equations, for quasi-linear differential (and time-strain separable single integral) equations, the local stability analysis with the method of frozen coefficients can be employed without loss of generality [83]. The studies of Hadamard stability have also a long history in the theory of nonlinear elasticity, where many general results were obtained and understood in great detail. Among the many general stability conditions suggested, the simplest is known as the Baker-Ericksen inequality [84]. Its physical sense is that the maximal principal stress always occurs in the direction of the maximal principal strain. In the case of hyperelastic solids, the thermodynamic stability criteria called "GCN + conditions" were also established long ago (see, e.g., Section 52 in Ref. [4]). Their physical sense is the convexity of the elastic potential with respect to the Hencky strain measure. This has been also known as a condition for strong monotonicity of stress with respect to strain. In nonlinear elasticity, the Hadamard stability criteria of field equations correspond to the conditions of strong ellipticity which coincide with the stability requirements known for dynamic problems. The weaker conditions of ordinary ellipticity deliver the stability constraints for static states [85].

542 The GCN + condition is closely associated with the condition of strong ellipticity. The lack of symmetry in the representation of second-rank variables causes a more restrictive condition for the strong ellipticity than for the GCN + condition [4]. Therefore, such inequalities as the Baker-Ericksen and the GCN + may be treated as necessary conditions for the strong ellipticity or the Hadamard stability. For the important case of isotropic incompressible hyperelastic solids, Zee and Sternberg [85] have recently found the necessary and sufficient conditions for the strong ellipticity (Hadamard stability) of CEs in a form of algebraic inequalities imposed on the first and second derivatives of the elastic potential. For the compressible case, these results have been obtained earlier [86]. For viscoelastic CEs, general results on global Hadamard stability, i.e. stability for any type of flow and for any Deborah number, were recently obtained for both general classes of quasilinear Maxwell-like CEs [36] and time-strain separable single integral CEs [87]. It was found later that these algebraic criteria for Hadamard stability were confused with the necessary and sufficient conditions for thermodynamic stability, that is, GCN +. These convexity conditions for thermodynamic potential in the hyperelastic case impose weaker constraints on CEs than the criteria for Hadamard stability. The complete results on the Hadamard instabilities for both the classes of viscoelastic CEs under study were then found and published in paper [88] for incompressible viscoelastic CEs. For compressible case, the complete criteria for Hadamard stability of viscoelastic CEs were recently found by Kwon [89]. In contrast to the Hadamard instability, another important type of instability, of dissipative type, results from the formulations of the non-equilibrium (dissipative) terms in CEs. This type of instability can occur even if the dissipation is positive definite. The studies of dissipative instabilities, initiated in paper [36] for general Maxwell-like CEs, were motivated by the fact that the upper convected Maxwell model which is globally Hadamard stable, displays the unbounded growth of stress in simple extension when the elongation rate exceeds the half of the reciprocal relaxation time. For any regular flow with a given history, a sufficient condition for the dissipative stability, close to the necessary one, was proposed in [36] for CEs of the differential type. Then the necessary and sufficient condition for single integral CEs was found in paper [90]. It was also noticed in those papers that in many viscoelastic flows, neither strain nor stress histories are known, but rather a complex mixture of both. Several patterns of pathological behavior as predicted by some popular specifications of Maxwell-like CEs, related to dissipative instability, were exposed in 1D flows [91 ]. It was found that it is necessary for the dissipative

543

stability that both the steady flow curves in simple shear and in simple elongation have to be monotonously and unboundedly increasing. The combined Hadamard and dissipative stabilities have also been analyzed in paper [88] for incompressible viscoelastic CEs. The results of this analysis will be discussed at the end of this Section.

3.2 Hadamard stability criteria for viscoelastic constitutive equations We now briefly outline the mathematical procedure [88] to obtain the conditions for global Hadamard stability for quasi-linear differential and timestrain separable single integral CEs in incompressible isothermal flow of viscoelastic liquid. It has been proved several times [36,75,79,80] that the evolution equation (19) with the mixed convected time derivatives (~=+1) is Hadamard unstable. Therefore, we consider below only the case of upper convected time derivative in the evolution equation (19) of the differential type and single integral CEs of the form (36). Then the total set of equations in this stability problem consists of viscoelastic CEs (48)-(50) complemented by the equations of momentum balance and incompressibility: pdv_/dt =-Vp +V-o,

V.v=O.

(51)

We assume that the set has a solution {__c,v,p} which satisfies some proper initial and boundary conditions, and impose on the solution extremely short and high frequency, infinitesimal waves of disturbances" {8g, 8g, 8v, 8p} = ~{ ~, g, v, p }.exp[i(k.x-mt)/g 21

(52)

Here c is a small amplitude parameter (leJ << 1), c, v, p and p are (generally complex) amplitudes of the corresponding disturbances, k is a wave vector, and co is the frequency. Using just the local linear stability analysis, we can easily find the following "dispersion relation", i.e. the dependence of the frequency co on the wave vector k and the parameters of the basic flow: for differential models

1 ~,)2~2 _ I B~j'n~V~kjVmk~ 2

[_~m(t - t 1)Bij,~mdt 1 9v i k j V m k n kjvj -0.

for integral models

(53)

544

Here ~ = co- k.v is the frequency of oscillations, with Doppler's shift on the basic velocity field v taken into account. The last equation in (53) means that due to the incompressibility condition in (51), the wave vector of disturbances k is orthogonal to the amplitude of velocity dis.turbance v. The fourth rank tensor Bijmn is defined in the principal axes of tensor __cwith principal values ck as: Bijnm- (~,n~n + ~ln~,n )Gij + ~j~nnLim, Gij = (Ci q- Cj)[(q)l q- q)2(I1 - C i - Cj)],

(no sum!)

(54)

Lij = 4q)nfifj + 2q~12CiCj(211 - Ci- Cj) - 2q)21(CiCj-1+ ci'lfj) Formulae (54) show that the tensor Bjjmn depends linearly on the constitutive functions q0k and their second derivatives q0kj with respect to the basic invariants Ik, independently of whether the approach is potential or not. It is now seen that due to equation (52) and the definition of ~ , the requirement for the stability is that the left-hand side of equation (53) should be positive. Therefore, the necessary and sufficient condition for global Hadamard stability is reduced to the following: 1 ~-~2V2

2

- BijnmvikjVmk n > 0

(55)

For the integral type constitutive equations, the condition (55) provides the stability during relaxation, which is included in the global stability requirement. On the other hand, the convexity of potential F, or the GCN + condition, can be represented as follows"

B ~jmn[30]3,,m > 0,

-02 F B ijmn -- Oh ~j0h,m

= 4Cmq

0c qn

eip

0C pj

Here hij is the Hencky strain measure, h__ = (1/2) I n c . To guarantee the thermodynamic stability, the inequality in (56) should be satisfied for any arbitrary symmetric tensor ,fl,j with the condition of incompressibility, trl~ = 0. In the potential (hyper-viscoelastic) case, the identity, Bijnm = B ijmn, holds. The comparison between the inequalities (55) and (56) shows that due to the symmetry of the tensor fl,j, the condition (56) is included in the inequality (55), i.e. the condition (56) imposes weaker stability constraints than the inequality

545

(55). It means that the conditions for the Hadamard stability are stronger than those of GCN +. Employing the algebraic procedure which has been used in hyper-elasticity [85], one can finally obtain the necessary and sufficient conditions for the global Hadamard stability as the set of algebraic constraints imposed on the functions q~k: (i) lai > O,

~i = ((Pl + (p2Ci)~k/CjCk

(i~j~k),

(ii) r q- 21-ti> O, q i =(Ii-Ci)(qOl+q02Ci)+2(Ii2-212-Ci 2-

(iii) [

,+21ai +

21

q>O,

s. )[q)llW(q)12+q)21)ci-l-q)22Ci 2] Ci

(57)

j + 2 ~ j ] 2 > g k - 2 ~ t k (iCj~k).

The additional constraint in (i), the positive definiteness of the tensor c, holds by definition for the integral CEs. It was also proved for the Maxwell-like CEs of differential type [36,102], however, only for the flow situations with a given history. The above approach to the global Hadamard stability has been recently extended by Kwon [93] on the compressible case. The new quality which occurs there is the possibility of longitudinal wave propagation. In the incompressible case, the speed of the longitudinal wave approaches infinity, whereas the speed of the transverse wave is finite. Hence, perturbation of basic solutions by the longitudinal wave was not considered in this stability analysis. The result was that the wave vector is always orthogonal to the vector of the main velocity field. However, in the compressible case, the speeds of both waves have finite values. Thus for stability, the initially infinitesimal amplitude of disturbing waves of either type (or a mixed type) should remain small all the time. It means that the conditions of Hadamard stability in compressible case are more rigid than those for incompressible one. It was demonstrated [89] on the simple example of Mooney-Rivlin potential with additional term dependent on density. Two sufficient conditions for the incompressible case have also been proposed: (1) The author's condition (27): the thermodynamic potential F for the author's class of viscoelastic CEs is a monotonously increasing convex function of invariants I1, and I2.

546 (2) Renardy's condition" the thermodynamic potential F for the K-BKZ class of CEs is a monotonously increasing convex function of ~ and Although Renardy's condition has been proved only for the K-BKZ class, it also holds for the Maxwell-like CEs with upper convected time derivatives [88]. Since the author's condition (1) is stronger than Renardy's (2), it also guarantees the global Hadamard stability for the K-BKZ class. The above sufficient conditions are much more easier to employ than the necessary and sufficient conditions for Hadamard stability (57). Therefore they are very useful for a brief evaluation of the stability for new formulations of CEs. For the compressible CEs, one sufficient condition for the global stability has also been suggested [89], but it is too complicated to use.

3.3 Dissipative stability criteria for viscoelastic constitutive equations As mentioned, there can be another source of instability originated from specification of dissipative terms in viscoelastic CEs. For viscoelastic CEs of differential type, this instability may happen due to an improper formulation of the dissipative term ~ (or ~b when ~ = 1) in equations (19), even for the Hadamard stable CEs with positively definite dissipation. For single integral CEs, the instability results from fading memory effects in equations (30) and (36). Although the global criteria for dissipative stability of viscoelastic CEs are far from being complete (if it is in general possible), we discuss in this Section two specific criteria that have been proven. In the case of compressible flow, no theorem on dissipative stability is known yet, but the following theorems are presumably valid also for the compressible CEs.

3.3.1 Criterion I of dissipative stability Theorem 1.1 (the case of CEs of the differential type [36]). Consider the set of upper convected Maxwell-like CEs (8) with the positive dissipation D = D(T, Ii, 12, 13) defined in equation (9). Let the free energy F be a non-decreasing smooth function of three invariants Ik. If for any positive number E, the asymptotic inequality

O > E'll ell

when Ilcll oo

(11 11- (trc2)m)

(58)

holds, then in any regular flow, the configuration tensor =e and the stress tensor ere are limited.

547

Theorem 1.2 (the case of single integral CEs [90]). In any regular flow, the functionals of free energy (35) and dissipation (37) are bounded, if (and only if) the thermodynamically or Hadamard stable potential function F(H1,H2,H3), expressed in terms of principal Hencky strains Hk, increases more slowly than exponentially. In the theorem 1.2, principal values of Hencky strain tensor and Finger tensor for the total deformation are related as" Hi =

(1/2)lnCi,

or

=/1= (1/2)lnC,

trH = 0.

(59)

Detailed proofs and definitions are given in the papers [36,90]. While the theorem 1.1 has been proved for differential CEs as a sufficient condition close to the necessary one, theorem 1.2 provides the necessary and sufficient condition for boundedness of single integral CEs. The above theorems were motivated by the fact that the globally Hadamard stable upper convected Maxwell model displays the unbounded growth of stress in simple extension when the elongation rate exceeds the half of the reciprocal relaxation time. As the consequences of the above theorems, (i) the upper convected Maxwell model which violates Criterion I, and (ii) the K-BKZ class with a potential F represented as an increasing rational polynomial function of basic invariants Iu, are dissipative unstable. Therefore, the Mooney and the neoHookean potentials as well as the potentials for the K-BKZ class of CEs which are subordinate to Renardy's sufficient evolution criterion also violate Criterion I of dissipative stability. Since the satisfaction of only Criterion I cannot prevent the viscoelastic CEs from severe dissipative instability, an additional criterion for dissipative stability has been introduced.

3.3.2 Criterion II of dissipative stability [88] For the stability of Maxwell-like and time-strain separable single integral CEs, it is necessary that both the steady flow curves in simple shear and in simple elongation have to be monotonously and unboundedly increasing with respect to the strain rate. It has been demonstrated [91] that the violation of Criterion II results in "blow-up" instability or even negative principal values of tensor __c in simple shear. Therefore the subordination to the combined criterion "I+II" was assumed in [88] to be presumably sufficient for the dissipative stability of both the differential Maxwell-like and the time-strain separable single integral CEs, at least in simple flows.

548 3.4. Application to viscoelastic CEs. Discussion Both the Hadamard and dissipative types of instability for such two broad classes of viscoelastic CEs have been discussed in this Section. These are the quasi-linear differential and factorable single integral models with instantaneous elasticity, which are the only ones in practical use today. The problem of global Hadamard stability for these two classes of CEs seems to be completely resolved in the isothermal, incompressible and compressible cases. This problem was reduced to that well known in the nonlinear elasticity, where the complete set of necessary and sufficient conditions of stability was formulated in algebraic form. It has been demonstrated that the proposed analysis of Hadamard stability for the two classes of viscoelastic liquids is reduced to the analysis of strong ellipticity. The physical sense of this is very evident: the studies of Hadamard stability involve very rapid disturbances which create only elastic response in viscoelastic liquids. In the case of dissipative stability, the global analysis is far from being completed, if it is generally possible. However, two distinct patterns of dissipative instability have been revealed, which are related to (i) the boundedness of stress, free energy and dissipation in a start-up flow problem under a given strain history (Criterion I), and (ii) the monotonously and unboundedly increasing steady flow curves in simple shear and simple elongation (Criterion II). Furthermore, it was assumed that the subordination of CEs to the combined criterion "I + II" is presumably sufficient for the dissipative stability in the simple flows. There is a tough problem as to how to distinguish the unstable behaviors caused by poor modeling of CEs and the observed physical instabilities which those equations should also describe. However, the long history of various branches of continuum mechanics and physics teaches us that the occurrence of either Hadamard instability or/and ill-posedness in ID situations without such important physical reasons as phase transitions, etc., is a distinct sign of inappropriateness in the CEs. Thus we can treat the instabilities demonstrated in this section as being associated not with the real instabilities observed in flows of polymer melts, but rather with the improper modeling of various terms in CEs. In numerical simulations of complex flows with unstable CEs, when the flow rate becomes high enough, the occurrence of various types of unphysical instabilities is inevitable. Even in the range of moderate Deborah numbers, the existence of singular points in flow geometry such as the comer singularity in die entrance region, is sufficient to spoil the entire numerical procedure.

549 All the results of the stability analyses found in various studies for popular viscoelastic CEs, are summarized in Table 1. An interested reader can find the details of calculations in references also provided in the Table 1. It is noteworthy that CEs derived from molecular approaches such as the Larson and the Currie models, exhibit the most unstable behavior. Surprisingly enough, none (to the authors' knowledge) of the time-strain separable single integral models are evolutionary. Appendix A represents the explanation of the reasons for that given by Simhambhatla [94]. He analyzed the time-strain separability concept for CEs and concluded that the Hadamard unstable CEs of time-strain separable type cannot properly describe the experimental data of stress relaxation after step-wise loading. The instabilities revealed in Ref. [94] exactly correspond to the results reviewed in this Section. It is astonishing that many CEs become Hadamard unstable even in viscometric flows. For the CEs of differential type, only three stable specifications exist. These are the FENE, the upper convected Phan-Thien-Tanner models, and the author class of CE's (8), (26) under convexity constraints (27). However, the FENE and the upper convected Phan-Thien-Tanner models predict zero value for the second normal stress difference in simple shear flow, which contradicts the experimental evidence for polymer melts and concentrated polymer solutions. It should be noted that all the necessary and sufficient conditions obtained for single-mode CEs become, strictly speaking, only sufficient for the multi-modal approach. Even though the necessity is not proved, it is thought that due to the inequalities (29), i.e. well separateness in the relaxation times, the exact conditions for Hadamard stability exposed above for a single mode CE, will be closed to necessary for multi-mode approach. It is also evident that the threshold of instability would only be delayed to some higher Deborah number region in the multi-modal approach, if any single-mode is unstable. For some viscoelastic CEs, regularization of ill-posedness may be achieved. E.g., it is well-known that adding a small Newtonian term to the stress stabilizes Hadamard unstable CEs. However, for complex flow simulations, this may not be enough to suppress numerical instability, and when the Newtonian term becomes larger, the description of the CE will deviate from the experimental data. In the case of Hadamard stable but dissipative unstable CEs which violate the Criterion II, one can also propose the more fundamental procedure of stabilization by changing the elastic potential. For example, the Giesekus model with or
550 Table 1 shows that the combined stability criteria impose very tough constraints on viscoelastic CEs. Therefore, the serious question arises as to whether there exists a CE or a class of CEs which can properly describe all the available rheometric data for concentrated polymer solutions and melts, when satisfying all the stability constraints. We will describe such a class in the next Section. Table I Stability of viscoelastic constitutive equations Type of instability Model (Eq. #) Type of CE Dissipative Upper convected Quasilinear unstable(Criterion I) Maxwell (20) (~=1) differential Hadamard unstable Interpolated Quasilinear Maxwell(20) differential (~,1) (Johnson/Segalman Gordon/Schowalter) Hadamard unstable General Phan-Thien/ Quasilinear (~.l) Tanner (21) differential Hadamard stable; Upper convected Quasilinear dissipative stability Phan-Thien-Tanner differential depends on (21) (~=1) dissipative term Nonlinear differential Hadamard unstable, White-Metzner(22) dissipative unstable (Criterion I) Globally Hadamard FENE(23) Quasilinear and dissipative stable differential Dissipative unstable Giesekus(24) Quasilinear (Criterion II) differential Dissipative unstable Simplest Leonov(26) Quasilinear (Criterion II) (bl, = b 2 = 1) differential Globally Hadamard Leonov class (26) Quasilinear dissipative stable under stability and differential constraints (27) Hadamard unstable Larson (25) Quasilinear differential dissipative unstable (Criterion II) single Hadamard unstable Wagner I (38) Separable integral single Hadamard unstable Wagner II (39) Separable integral single Hadamard unstable Papanastasiou (40) Separable (K=0) integral

References e.g.[36] [36,75,76,80]

[36,75,76,80] [88]

[75,81] [88] [91] [91] [88] [88] [91] [88] [88] [90]

551 Luo-Tanner (40) Lodge(41) K-BKZclassunder Renardy's condition Larson-Monroe potential (42)

Separable single integral Separable single integral K - B K Z Separable K-BKZ Separable K-BKZ

Currie potential (43)

Separable K-BKZ

Yen-Mclntire (44)

Quasi-separable K-BKZ

Hadamardunstable

[88]

Dissipative unstable [90] (CriterionI) Dissipative unstable [90] (Criterion 1) Hadamardunstable [88] dissipative unstable (Criterion II) Hadamard unstable dissipative unstable (Criterion II) Dissipative unstable (Criterion 1)

[90] [88] [95] [90]

4. MODELING OF POLYMER FLUIDS WITH STABLE CONSTITUTIVE EQUATIONS Following mostly the paper [40], we demonstrate in this Section the class of differential CEs [23,24,36] which is able to consistently describe simple flow data for such basic polymers as HDPE II, PS I, PIB P-20, PIB L-80 and LDPE Melt I/ IUPAC A/ IUPAC X, while complying with the global isothermal stability constraints. For simplicity, only one or two nonlinear parameters additional to the discretized linear viscoelastic spectra are introduced for the description of data. Instead of the simple "Leonov model" which uses only the parameters of the discretized linear viscoelastic spectrum, we employ in the following sections a highly nonlinear specification of the general class of Maxwell-type CEs proposed by the author [23,24,36]. This specification subordinated to the convexity conditions (27) guaranties the both Hadamard and dissipative stabilities. Comparison of the descriptive ability of other models with experimental data is not attempted in this Section. 4.1. Selection o f a descriptive subclass from the author's CEs For the sake of simplicity, only the constitutive modeling of viscoelastic liquids with incompressible Maxwell-type CEs is considered below, with the equations shown for a single relaxation mode. The model equations are given by equations (2), (8), (26) and (27).

552

We consider firstly the modeling of dissipative terms in evolution equations (8) and (26). The following simple forms of equation (26) have been considered

[40]: i) b~ = b(I~,/2,T)/20(T), b2 = 0; ii) b~ = 0,

b2 = b(I~,/z,T)/20(T);

(60)

iii) b~ = b2 = b(I~,Iz,T)/20(T). Here, 0 = 0(T) is the relaxation time in the linear Maxwell limit. The specification i) with the Neo-Hookean elastic potential, results in a decreasing branch of the flow curve in simple shear for b = 1. This poor quality, resulting in dissipative instabilities, can be rectified by specifying more sophisticated functional dependencies either for b(I~,I2,T), or the elastic potential. However, this form was rejected due to the inconvenience for practical modeling. With the specification ii), various simple flows can be described quite accurately for several polymers. The only problem with ii) is the weak maxima of Nl predicted during start up shear flow, in comparison with experimental observations. No way of rectifying this malady was found and we therefore reject this specification also. The form iii) for which the evolution equation (8), (26) can be written as: v

20(T)c +b(I, ,12,T) [c2 + c.(I2-I,)/3 - __6]= __0

(61)

is the only one discussed and tested in the literature. It allows for plane deformations in simple shear and endows the resulting equations with the proper quality and flexibility for modeling a wide variety of polymers. Detailed comparisons with experimental data have been made [40] with this representation. In order to relate the extra stress tensor to the elastic Finger tensor during the deformation history, a functional form for the elastic potential W(I1,/z,T) = poF must be provided. Here, Po is the density, and F is the specific Helmholtz free energy. The following fairly general elastic potential:

3G(T)

W(II,I2,T) - 2(n+ 1) {( 1-[3)[(I~/3)n+1-1] + (1-13)[(I,/3)n +'-1 ],

(62)

553

has been suggested in [40]. Here G(T) is linear Hookean elastic modulus, and 13 and n are numerical parameters. With the evolution equation (61) and potential relation (62), the constitutive equations are Hadamard and dissipative stable for 0_< 13 <1 and n >0 (see the sufficient condition O for Hadamard stability and criterion 'I +II' for dissipative stability) except when the criteria for "fluidity loss" are met (see below). Equation (62) yields the Mooney potential for n = 0, and the neo-Hookean potential for n = 13= 0. The extra stress tensor can then be written due to equation (18) as: flex = (1-13) (I1/3) n __C-l] (12/3) n C"1.

(63)

In the multi-mode approach employed below the same values of nonlinear constants for each mode are used. This gives a few-parametric description of the data. Also, the conditions of Hadamard and dissipative stability used for every mode will obviously satisfy the stability criteria for the complex model. For modeling isothermal experimental data, we need specifications for the function b(1t,12), as well as a simple reduction of the potential W(It, I2) suggested in equation (62). For practical modeling, b(I1,/2) can be thought of as a deformation-history dependent scaling factor for the linear relaxation times. The simplest choice is to let b = 1, which is known as the standard "Leonov model". Henceforth, the CE with b=l and the neo-Hookean potential (n = 13= 0) is further referred as the "simple model". While this simple choice assures the proper linear viscoelastic limit, and can also be expected to describe weak nonlinearities, it may not suffice for the description of highly nonlinear phenomena. For instance, some polymers (e.g., LDPE with a high degree of branching) show great hardening in simple elongation, while others (e.g., HDPE, PS) do not. In any simple flow, if b(Ii,/2) is chosen to decrease gradually (e.g., using weak power law) with an increase in the magnitude of 11 and 12, there will be hardening relative to the simple model. A rapid fall in b(Ii,/2) (e.g., exponential) will cause the steady state components of c to be double-valued up to a critical value of the strain rate, with one stable branch. Beyond this critical value, there is no steady-state solution and the components of c increase unboundedly as in an elastic solid. This is the concept of "fluidity loss" analyzed in detail by the author [23,32]. On the other hand, an increasing b(Ii,/2) will cause relative softening. The more rapidly b increases, the more gradual will be the variation of the steady-state value of c_with an increase in the strain rate. The choice of b(I~,/2) is not as difficult as it appears. The recommended procedure is to first perform some preliminary calculations for various flows

554

with the simple model. Then, if there is disagreement with experimental data, an appropriate functional form of b function required to bring the calculations into qualitative agreement with the data, can be systematically determined. The reason this procedure is straightforward is that for this class of equations, the effects in various flows are well separated in the sense that there is some measure of flexibility in modeling their qualitative behavior independently. This is a direct consequence of the fact that whether the simple deformations are steady or non-steady, the following relationships hold for the invariants of the elastic Finger tensor c: I~ = 12

in simple shear and planar elongation

11 > 12 in simple elongation, and

(64)

11 < I2 in equi-biaxial extension.

These relations also hold for the invariants of the total Finger deformation tensor C. This is a remarkable feature of the evolution equation (61). The experience [40] in modeling the viscoelastic behavior of several polymers (LDPE, HDPE, PS, low and high molecular weight PIB) with the function b(I~,/2) showed that simple power law or exponential functions of the invariants with one or two adjustable parameters are sufficient for an accurate quantitative description of all the available data. However, unlike the elastic potential (62), no any unified form for b(I1,/2) has been found which could be used for the description of data for all polymer melts. Because of the flexibility which the modeling of the dissipative term permits, one can operate with fairly simple forms of the general elastic potential (62), such as the neo-Hookean and Mooney potentials for the description of the usual rheometrical data. However, in simple shear, with this approach, cl2 generally reaches a limit value of unity at high Deborah numbers leading to a saturation of the shear stress and therefore to dissipative instability. One remedy is to extend the discrete relaxation time spectrum at the small relaxation time end in order to effect stability until the region where physical instabilities appear. A simpler approach is to work with the existing discrete relaxation modes obtained from the usual linear experiments, and to allow the parameter n in the potential relation (62), to be a small positive number (e.g., n = 0.1). In this way, we can preserve all the predictions in the rheometric regimes with simple potentials and also effect the non-saturation of the stress at very high Deborah numbers. Still, for concrete recommendations for quantitatively modeling high

555

Deborah number flows, more data in the region of incipient physical instabilities would be welcome.

4.2. Component equationsfor simple flows For convenience, the equations and initial conditions for simple flows for a single Maxwell mode are presented below. These equations and formulae should be employed in a multi-mode approach for all the predictive calculations in comparison with data.

4.2.1. Simple shear The evolution equations take the form" 20dCll /dt + b(I)(c~ + c 212- 1) - 4,~0c 12 (65)

20dc12 / dt + b(I)Cl2 (Cll +

C22 ) --

2~C22

Cl.C.., ..-l+c

-I-1

+ c 1 + C 1 22 9

212 9~

I 1 -I

2

The system of stresses is: cy12 - G(I / 3) n C12 ;

N 1 -

G(I / 3) n ( e l l

_ C2 2 ) ;

(66)

N 2 - G(I/3)n[(1 -- ]3)(C22 -- 1)+ ]3(1-- C11)]

Here, ~, is the shear rate, O"12 is the shear stress, Nl and N2 are the first and second normal stress differences, respectively.

Startupflow The initial conditions for startup flow are: cijl t=o- 8ij.

(67)

The steady-state solution of equations (65) is of the form" Cll -- - ~ Z / 4 Z + 1;

r

Z ( I ) - (I / 2)(I - 1) 2 - 1.

- ~

/ 4 z + 1;

Ol). - ~ / z ~ - 1 / ( z + 1); (68)

556 Here, I can be obtained by solving the implicit equation" Z ( I ) - 41 + (20~ / b(I)) 2 .

(69)

Stress relaxation Here, ~, = 0. For stress relaxation following cessation of shear flow, the initial conditions for equations (65) are:

c,j[ t_-o- cij]v,t ,

(70)

where ts is the shearing time prior to cessation of flow. For the stress relaxation following the imposition of a step strain 7o, the initial conditions for equations (65), are: ell

]t = O * - 1 + 7 2. o,

]

c22 t-o*

-1;

c 12 ]t=o* -70

(71)

9

Creep and recovery For creep, let a shear stress ,:3-~ be applied at time t = 0. Let the shear strain at time t = +0 after the jump be %. Then the initial conditions for the kinematic variables are given by equation (71). The value of ~,o is obtained by substituting for C12 from equation (71) into equation (66) and setting O"12 - - CY~ For the multi-modal case, ~ ~2 is the total stress defined as the sum of the sub-stresses in the various modes. Then the shear rate, ~; is calculated by setting dcrlz/dt = 0 in equation (66) and substituting for the time derivatives of the kinetic variables from equation (65). In the multi-modal case, r~2is the total shear stress. Evidently, d7 / dt - ~

(7] t=0" - 70).

(72)

For N modes, the set of 2N+ 1 differential equations, two per mode in equation (65) and one in equation (72), can then be solved by using e.g., a Runge-Kutta scheme with an automatically adapting step size. The conditions for recovery following unloading can be similarly derived.

4.2.2. Simple elongation The evolution equation is:

X-~dX/dt+b(I~,I2)(X 2-X-~)/(6~.0)-~;;

I~-~2+2X-~

12

_ ~-2 +2~.

557 (73) where ~ is the elongation rate. The elongation stress is: G"E

=

G [ ( 1 - 13)(1, / 3) n (X 2 - X-') + 13(I2 / 3)(X - X-2 )].

(74)

Startup flow The initial condition for startup flow is: Lit_ o - 1.

(75)

The steady-state solution for X is obtained from the implicit equation: (76)

b(I, ,I 2 )(~2 _ X-' )(1 - k-') = 6/~0.

Stress relaxation Here, ~=0. For stress relaxation following cessation of elongation flow, the initial conditions for equations (73) is:

~l t=0-- ~l+,t

,

(77)

where te is the time of extension before cessation of flow. For stress relaxation following the imposition of a step Hencky strain eo, the initial condition is: )~l t=o" - exp(eo)"

(78)

The procedure for creep and recovery calculations is analogous to that for simple shear.

4.2.3. Planar extension The evolution equation is: ~-' d~ / dt + b(I)(X 2 - ~-2) / (40(T)) = ~p;

I, - 12 - I - )~2 + )C2 + 1.

(79)

Here ~p is the planar extension rate. The system of planar extension stresses is

558

o p, - 0,, - ~22 - G(I / 3)" ( 1 - 213)(X2 - )v-2) 0 p, - 0 22 - 0 33 - G(I / 3)" [13)v2 + ( 1 - 13))v2 - 1].

(8o)

Startup flow The initial condition for startup flow is" )q t - o - 1 .

(81)

The steady-state solution for 9~ is obtained from the equation: b(I)(L 2 - X-2) - 4~p0.

(82)

Stress relaxation Here, /~p= 0. For stress relaxation following cessation of planar extension flow, the initial condition for equation (79) is:

~1t=0-- ~l tp,~p,

(83)

where tB is duration of biaxial extension before cessation of flow. The initial condition for stress relaxation following the imposition of a step biaxial Henckey strain eBo may be written as: ~1 t=O -- exp(e Po ).

(84)

Creep and recovery calculation can be performed as for the other simple flows.

4.2.4. Equi-biaxial extension The evolution equation is: L-'dX / dt + b(I,,I2)(X 2 - ~-,)(~2 + 1) / ( 1 2 0 ) - ~B; I~ - 2X 2 + X -~"

I2 - 2X-2 + X~

~

(85)

Here, ~ B is the biaxial extension rate. The biaxial stress is" ~B - G[(1 - 13)(1, / 3) n (X: -- ~-4 ) ..1_ ]~(i 2 / 3) n (~4 __ X-2 )].

(86)

559 Startup flow The initial condition for startup flow is"

)q t=o- 1.

(87)

The steady-state solution for L is obtained from the equation: b(I,,I2)(~, 2 - X-4)(X2 + 1 ) - 12/~,0.

(88)

Stress relaxation Here, ~B= 0. For stress relaxation following cessation of biaxial flow, the initial condition for equation (85) is"

~] t=0-- )L]g,,tB,

(89)

where tB is duration of biaxial extension before cessation of flow. The initial condition for stress relaxation following the imposition of a step biaxial Henckey strain eBo may be written as" XIt-o - exp(~:Bo)"

(90)

Creep and recovery calculation can be performed as for the other simple flows.

4.3. On the comparison with experimental data High density polyethylene HDPE-II, polystyrene PS-I, polyisobutylene (PIB) P-20, a relatively high molecular weight PIB, Exxon Vistanex L-80, and low density polyethylene LDPE Melt-I have been chosen in paper [40] to compare the predictions of the above constitutive equations with data. A numerous amount of data and calculations involved in the comparison demonstrated generally a great success in our modeling. The interested reader can find a lot of useful details in Ref.[40]. For four first polymers in the tested group, the specification of the function b(Ii,I2) was uniform and proposed as follows: b(h) = exp[m(/j3 - 1)].

(91)

560 Eq.(91) means that the four first tested polymers demonstrate softening behavior at high strain rates. For describing the data for LDPE Melt I with the author's class of CEs, some preliminary calculations were first performed with the simple model, as suggested in the section 3.1 [40]. With this approach, the description of shear flows was quite accurate. However, the biaxial extension damping function was overpredicted, while the hardening effects in extension flow were underpredicted. To rectify the observed discrepancies with the simple model, the relations (60) suggest that b(11,I2) should be a decreasing function of (IJlz) (see section 3.1 for the physical sense of the function b). A simple choice was made as: b(11,I2) = (I2/11) '~ .

(92)

The parameter 'm' was chosen to be 1.4, for properly describing the extension stress growth data. It was a hope [40] that with formula (93) for b(I1,I2) all the available experimental data could be described reasonably well. Indeed, the calculations according to this choice of b(I1,I2) could describe properly almost all the data for the Melt I but they failed to describe the planar elongation tests [107]. The reason for this was that in the planar elongation, as in the simple shear, 11 = 12. Thus the hypothetical rheological behavior in the planar elongation, as predicted by formula (92), is softening. However, this prediction contradicts the hardening phenomena in planar elongation, observed experimentally [ 107]. Therefore in Appendix B, a new, more physically related attempt is presented to describe the whole set of data for LDPE Melt I.

5. CONCLUSIONS The behavior of two common classes of viscoelastic constitutive equations (CEs) for polymer melts and concentrated polymer solutions was discussed. These are general Maxwell-like and single integral CEs with instantaneous elasticity. The formulation of both classes of CEs was analyzed. The Maxwell-like CEs usually employ some hidden tensor variables with different physical senses. Therefore, in spite of the generality in formulation, their evolution equations and stress relations have different features, depending on the theoretical approach used. Some artifacts related to formulation of the CEs were also exposed. Such an important effect as compressibility was discussed.

561

General results on stability for both classes of CEs were demonstrated, which included stability analyses of both the Hadamard and dissipative types. Results of the stability analyses were applied to popular CEs. The descriptive capability of a class of Maxwell-like CEs was demonstrated, whose formulation satisfies all the stability constraints. It should be noted that the data [98] for equi-biaxial extension were obtained with using lubricated squeezing technique. Recent publications [113,114] reported that experiments with this technique can involve undesirable side effects, such as distortion of sample shape [113] or uncontrolled thinning of lubricant layer [114]. Therefore there is still a need for independent equibiaxial extension data for polymers. REFERENCES ~

2. 3. o

.

6. 7. 8. ,

10. 11. 12. 13. 14. 15. 16. 17. 18.

Odroyd, Proc. Roy. Soc., A200 (1950) 1063. Oldroyd, Proc. Roy. Soc., A245 (1958) 278. Rivlin, Research Frontiers in Fluid Dynamics, 144, Interscience, New York, 1965. Truesdell and W. Noll, The Non-Linear Field Theories of Mechanics, Springer, New York ,1992. Coleman, Arch. Rat. Mech. Anal., 17 (1964) 1. Coleman, Arch. Rat. Mech. Anal., 17 (1964) 230. Coleman and V.J. Mizel, Arch. Rat. Mech. Anal., 29 (1968) 105. de Gennes, Scaling Concepts in Polymer Physics, Cornell University, Ithaka, New York, 1979. Green and A.V. Tobolsky, J. Chem. Phys., 14 (1946) 80. Lodge, Trans. Faraday Soc., 52 (1956) 120. Yamamoto, J. Phys. Soc. Japan, 11 (1956) 413; 12 (1957) 1148; 13 (1958) 1200. Treloar, Physics of Rubber Elasticity, 3rd ed., Oxford University Press, Oxford, 1975. Edwards, Proc. Phys. Soc., 92 (1967) 9. Doi and S.F. Edwards, J. Chem. Soc., Faraday Trans. II, 74 (1978) 560, 918, 1818; 75 (1978) 38. Curtiss and R.B. Bird, J. Chem. Phys., 74 (1981) 2016. Bird and J.R. Deaguiar, J. Non-Newton. Fluid Mech., 13 (1983) 149, 161. Curtiss, R.B. Bird and O. Hassager, Adv. Chem. Phys., 35 (1976) 31. Bird, R.C. Armstrong and O. Hassager, Dynamics of Polymeric Liquids, Vol. I & II, 2nd Ed., John Wiley and Sons, New York, 1987.

562 19. 20. 21. 22. 23. 24. 25. 26. 27. 28. 29. 30. 31. 32. 33. 34. 35. 36. 37. 38. 39. 40. 41. 42. 43. 44. 45. 46. 47. 48.

Pokrovsky and V.S. Volkov, Visokomolekulamie Soed. (Russian), A20 (1978) 2700. Volkov, Int. J. Polym. Mater., 9 (1982) 115. Schweizer, J. Chem. Phys., 91 (1989) 5802, 5822. Schweizer, J. Non-Cryst. Solids, 131-133 (1991) 643. Leonov, Rheol. Acta, 15 (1976) 85. Leonov, J. Non-Newton. Fluid Mech., 25 (1987) 1. Dashner and W.E. Vanarsdale, J. Non-Newton. Fluid Mech., 8 (1981) 59, 69. Larson, Constitutive Equations for Polymer Melts and Solutions, Butterworth, Boston, 1988. Hunter and M. Slemrod, Phys. Fluids, 26 (1983) 2345. Ramamurthy, J. Rheol., 30 (1986) 337. Tanner, Engineering Rheology, Clarendon Press, Oxford, 1992. Janeschitz-Kriegl, Polymer Melt Rheology and Flow Birefringence, Springer, New York, 1983. Pearson, Mechanics of Polymer Processing (Part I: Polymers, Mechanics and Rheology), Elsevier, New York, 1985. Leonov and A.N. Prokunin, Nonlinear Phenomena in Flows of Viscoelastic Polymer Fluids, Chapman & Hall, London, 1994. Grmela, Physica, D21 (1986) 179. Grmela, Phys. Lett., A 130 (1988) 81. Beris and B.J. Edwards, J. Rheol., 34 (1990) 55,503. Leonov, J. Non-Newton. Fluid Mech., 42 (1992) 323. Jongschaap, Physica, A94 (1978) 531. Jongschaap, K.H. de Haas and C.A.J. Damen, J. Rheol., 38 (1994) 769. Kwon and S. Shen, Rheol. Acta, 23 (1984) 217; 24 (1985) 175. Simhambhatla and A.I. Leonov, Rheol. Acta, 34 (1995) 259. Leonov, Polym. Internat., 36 (1995) 187. De Groot and P. Mazur, Non-Equilibrium Thermodynamics, NorthHolland, Amsterdam, 1962. Gyarmati, Non-Equilibrium Thermodynamics. Field Theory and Variational Principles, Springer, New York, 1970. Hulsen, J. Non-Newton Fluid Mech., 38 (1990) 93. Murnaghan, Amer. J. Math., 59 (1937) 235. Leonov, Ann. N.-Y. Acad. Sci., 410 (1983) 23. Stickforth, Rheol. Acta, 25 (1986) 447. Leonov, On the constitutive equations for non-isothermal bulk relaxation, accepted for publication in Macromolecules, 1997.

563

49. 50. 51. 52. 53. 54. 55. 56. 57. 58. 59. 60. 61. 62. 63. 64. 65. 66. 67. 68. 69. 70. 71. 72. 73. 74. 75. 76. 77. 78. 79. 80.

Gordon and W. R. Schowalter, Trans. Soc. Rheol., 16 (1972) 79. Johnson, Jr. and D. Segalman, J. Non-Newton. Fluid Mech., 2 (1977) 255. Phan Thien and R. I. Tanner, J. Non-Newton. Fluid Mech., 2 (1977) 353. Phan Thien, J. Rheol., 22 (1978) 259. Larson and K. Monroe, Rheol. Acta, 23 (1984) 10. White and A. B. Metzner, J. Appl Polym. Sci., 7 (1963) 1867. Giesekus, Rheol. Acta, 21 (1982) 366. Larson, J. Rheol., 28 (1984) 545. Simhambhatla and A.I. Leonov, Rheol. Acta, 32 (1993) 259. Rivlin and K.N. Sawyers, Ann. Rev. Fluid Mech., 8 (1971) 17. Kaye, College of Aeronautics, Cranford, U. K., Note No. 134 (1962). Bernstein, E. A. Kearsley and L. J. Zapas, Trans. Soc. Rheol., 7 (1963)391. Kwon and A.I. Leonov, Rheol. Acta, 33 (1994) 398. Wagner, T. Raible and J. Meissner, Rheol. Acta, 18 (1979) 427. Wagner and A. Demarmels, J. Rheol., 34 (1990) 943. Luo and R.I. Tanner, Int. J. Num. Meth. Eng., 25 (1988) 9. Papanastasiou, L. E. Scriven and C. W. Macosko, J. Rheol., 27 (1983) 387. Larson and K. Monroe, Rheol. Acta, 26 (1987) 208. Currie, in G. Astarita, G. Marrucci and L. Nicolais, "Rheology", Vol.1, Plenum, New York (1980). Yen and L.V. McIntire, Trans. Soc. Rheol., 16 (1972) 711. Zapas, J. Res. Natl. Bur. Std., 70A (1966) 525. Ogden, Proc. Roy. Soc., A326 (1972) 565. Mooney, J. Appl. Phys., 11 (1940) 582. Valanis and R.F. Landel, J. Appl. Phys., 38 (1967) 2997. Blatz, S.C. Sharda and N.W. Tschoegl, Trans Soc. Rheol., 18 (1974) 145. Crochet, A.R. Davies and K. Walters, Numerical Simulation of NonNewtonian Flow, Elsevier, Amsterdam, 1984. Dupret and J.M. Marchal, J. Non-Newton. Fluid Mech., 20 (1986) 143. Joseph, Fluid Mechanics of Viscoelastic Liquids, Springer, New York, 1990. Rutkevich, J. Appl. Math. Mech., 33 (1969) 30, 573; 34 (1970) 35. Godunov, Elements of Continuum Mechanics, Nauka, Moscow, 1978. Joseph, M. Renardy and J.C. Saut, Arch. Rat. Mech. Anal., 87 (1985) 213. Joseph and J.C. Saut, J. Non-Newton. Fluid Mech., 20 (1986) 117.

564 81. 82.

83. 84. 85. 86. 87. 88. 89. 90. 91. 92. 93. 94. 95.

96. 97. 98. 99. 100. 101. 102. 103.

104. 105. 106. 107. 108. 109. 110.

Verdier and D.D. Joseph, J. Non-Newton. Fluid Mech., 31 (1989) 325. Kreiss, Numerical Methods for Solving Time-Dependent Problems for Partial Differential Equations, Les presses de l'Universite de Montreal, Montreal, 1978. Strang, J. Diff. Eq., 2 (1966) 107. Baker and J.L. Ericksen, J. Wash. Acad. Sci., 44 (1954) 33. Zee and E. Sternberg, Arch. Rat. Mech. Anal., 83 (1983) 53. Knowles and E. Sternberg, Arch. Rat. Mech. Anal., 63 (1977) 321. Kwon and A.I. Leonov, J. Non-Newton. Fluid Mech., 47 (1993) 77. Kwon and A.I. Leonov, J. Non-Newton. Fluid Mech., 58 (1995) 25. Kwon, J. Non-Newton. Fluid Mech., 65 (1996) 151. Kwon and A. I. Leonov, Rheol. Acta, 33 (1994) 398. Kwon and A. I. Leonov, J. Rheol., 36 (1992) 1515. Hulsen, J. Non-Newton. Fluid Mech., 38 (1990) 93. Renardy, Arch. Rat. Mech. Anal., 88 (1985) 83. Simhambhatla, The Rheological Modeling of Simple Flows of Unfilled and Filled Polymers, Ph.D. Dissertation, the University of Akron, 1994. Kwon, Studies of Viscoelastic Constitutive Equations and Some Flow Effects for Concentrated Polymeric Fluids, Ph.D. Dissertation, The University of Akron, 1994. Einaga, K. Osaki, M.Kurata, S. Kimura, and M. Tamura, Polym. J., 2 (1971)550. Laun, Rheol. Acta, 17 (1978) 1. Khan, R.K. Prud'homme and R.G. Larson, Rheol. Acta, 26 (1987) 144. Soskey and H.H. Winter, J. Rheol., 28 (1984) 625. Takahashi, K. Taku, and T. Masuda, J. Soc. Rheol., Japan, 18 (1990)18. Laun, J. Rheol., 30 (1986) 459. Leonov and A.N. Prokunin, Rheol. Acta, 19 (1983) 137. Laun, Stress and recoverable strains of stretched polymer melts and their prediction by means of a single integral constitutive equation. In" Rheology, vol.2, Plenum Press, New York (1980). Munstedt and H.M. Laun, Rheol. Acta, 20 (1981) 211. Vinogradov and A.Ya. Malkin, J. Polym. Sci. A-2, 2 (1964) 2357; 4 (1966) 135. Leonov and A.N. Prokunin, Rheol. Acta, 19 (1980)393. Laun and H. Schuch, J. Rheol., 33 (1989) 119. Wagner, J. Non-Newton. Fluid Mech., 4 (1978) 39. Laun, Rheol. Acta, 21 (1982) 464. Wagner and H.M. Laun, Rheol. Acta, 17 (1978) 138.

565 111. Giacomin, R.S. Jeyaseelan, T. Samurkas and J.M. Dealy, J. Rheol., 37 (1993) 811. 112. Meissner, Trans. Soc. Rheol., 16 (1972) 405. 113. Takahashi, T. Isaki, T. Takigava, and T. Masuda, J. Rheol., 37 (1993) 827. 114. Kompani, D.C. Venerus, and B. Bernstein, Development and evaluation of lubricated squeezing flow technique. In: Proc. XIIth Int. Congr. on Rheology, A. Ait-Kadi, J.M. Dealy, D.F. James, and M.C. Williams, Eds. August 18-23, 1996, Quebec City, Canada, p.754.

Appendix A On the invalidity of strain-time separability at quick time scales (M. Simhambhatla [94]) Strain-time separability, i.e. the factorability of the material response to nonlinear step strains into time and strain dependent parts, has been widely used as a convenient basis for the specification of viscoelastic CEs. The proponents of this assumption claim justification based on experimental observations. However, it is shown here that the same experiments for polymer melts and solutions require that the principle of strain-time separability be violated at small times following the application of step shear strains, in order to guarantee the Hadamard stability. According to strain-time separability, the stress response, to a step shear strain applied from the rest state is: o(t,y) : yh(y)G(t),

(A 1)

where h(7) is the shear damping function, and G(t) is the linear (Maxwellian) relaxation function. This means that the stress relaxation curves (log(o(t)) vs. log(t)) for various applied step strains will be parallel. This observation of seemingly parallel stress relaxation curves has been reported for several polymer melts and solutions (e.g., see [96], [97]). Consider strain-time separable CEs with a perfect elastic limit. For these equations, the step stress cl in response to a step strain 7, will be N

a(7 ) - 7h(7 ) ~ G, i=l

(A2)

566

Here, N is the number of Maxwellian modes, and Gi the linear relaxation modulus for the 'i'th mode. Now consider the experimentally determined dependence Th vs. T for LDPE Melt I [98] shown in Figure A1. The maximum in this dependence implies that 6 vs. T should also have a maximum (equation (A2)). Such maxima for Th vs. T appear for all the experimental data that we have come across (e.g., see [96-100]). It is however, easy to see that the decreasing branch of the dependence 8 vs. T, is unstable in the Hadamard sense. In order to clarify this, consider a Cartesian coordinate system (x,y) with the 'x' axis parallel to the direction of the shear displacements 'u', and the 'y' axis normal to the shear planes. The equation of motion for simple shear can then be written as" p0v/0t = 3cr/o~

(A3)

With v - Ou/c~, and 3' = Ou/0y, equation (A3) can be rewritten in the elastic limit as; /~U/C~ 2= M(]t)c32u/o~ 2

(A4)

Here 'u' is the displacement in the 'x' direction, and M(T) - c36/0y. Let uo be a basic solution satisfying equation (A4). Let us impose a small disturbance on the basic solution, so that (A5)

u - Uo + ~ fi exp[i(ky-mt)/a21

where ~ is the small amplitude parameter, 'k' the wave number, and co the frequency of disturbance. Substituting for 'u' from formula (A5) into equation (A4) and taking into account only the lowest order terms in g, yields" (1)2 =

M(To)k2/p.

( M(y o ) - d r ~

03, Ir,, )

(A6)

If M(y o) > 0, m is real and the disturbances in equation (A5) will not grow with time. However, if M(To) < 0, we have: co - +i~/(]m(y o)lk

=/p.

(A7)

567

v _L--

2.0

r - - .........

1.5

J~

I

......

",~

...........

I

. . . . . . .

[

.........

1" ....

MeltI "

0

1.0

o5 [ o

-----

Experirnenfal h(7)

=

Poinfs

,57exp(--.31T)+.43exp(--.10-67) L ...........

5

i0

I

....

~ 5.. Shear

t, . . . . . . . . . . . . . . . . . . .

20

Stroin,

L

25

................. J

........................................

50

55

-y

Figure A 1. Experimentally observed maximum in the plot 7h(7)vs. 7 for LDPE Melt I [96]. When the imaginary part of o is negative, formula (A5) indicates that even infinitesimal high frequency disturbances will grow in amplitude rapidly and unboundedly in time. Because there are high frequency disturbances in the spectrum of natural noise, the decreasing branch of the dependence t~ vs. 7 will result in severe instabilities during experimental measurements. However, since the experiments have been conducted for large step shear strains, and no instabilities have been reported, it is evident that the principle of strain-time separability should be violated at small times following the imposition of these large shear strains. This violation should be manifested in the form of an upturn in the dependence cy(t) vs. t at small times following the imposition of a step strain 7 when moving along the time axis in the negative direction, in order to have a monotonic increase in the dependence ~ vs. the step strain Y. Interestingly, this is preciselywhat was observed by Einaga et al. [96] (see Figure 2) and Takahashi et al. [100]. In experiments where yh(7) vs. 7 has a maximum, but no loss of parallelism is seen in the curves ty(t) vs. t , it

568

is presumably because data could not been obtained at very short times following the imposition of the step strains (assuming no wall slip).

........

I ..........

~ ........

~-

10 4

t0 3 0,,,. v o p..

to2

v

r~ "10 (3)

,J~ r

101

i 011

.................................

1

I

1 0 `=

..............

. . . I ...... .. ................

10 2

i

10 s

.. . . . . . .

10 4

Figure A2. Illustration of violation of strain-time separability at short times for 20% PS in Arochlor [96]. It should be noted that the deviations from the master curve of vertically shifted stress relaxation curves are consistent with the requirement of Hadamard stability. Two conclusions can be drawn. The first is that strain-time separability holds only in an approximate sense, if we are prepared to neglect the material response to quick disturbances. The second is that any strain-time separable CE is certainly unstable in the Hadamard sense in the limit of very rapid disturbances, due to the universally observed maximum in the dependence ~,h(~/) vs. ~,in those time regimes. This indicates that we cannot simultaneously describe the experimental data in the rheometrical regimes, and have the stability at high Deborah numbers for any strain-time separable CE. It is evident that the above instabilities have no physical basis, but are simply rooted in the improper extrapolation of strain-time separability to quick time scales. High Deborah number flow simulations using these CEs should therefore be avoided.

569

Appendix B On rheological description of LDPE Melt I by stable CEs M. Simhambhatla and A.I. Leonov According to Exxon material data, LDPE Melt I has molecular weight Mw=460000 and polydispersity ratio Mw/Mn=22. Linear viscoelastic spectrum for this material at 150~ is shown in Table B 1. Table B 1 Linear__v_i.sc0e!ast_.!c_s_pectrum o f,L Dp E Melt i .(.15oOc)[9:7] i

1

2

3

0i, sec

10 3

102

101

G,, Pa

1.00xl0 ~

1.80x102

1.89x103

4 10 ~ 9.80x103

5

6

7

8

101

10 .2

10 .3

10-4

2.67x104

5.86x104

9.48x104

1.29X!0!

In order to describe all available data for the LDPE Melt I we will use the simple non-Hookean potential and two-parametric expression for the dissipative term b" b(I, ) - exp[-13(I, - 3)] +

sinh[v(I, - 3)1 vO, - 3 )

-1.

(0<13 < 1, v > O)

(B1)

Here [3 and v are some numerical fitting constants. Formula (B 1) is in fact a modification of a similar expression proposed in paper [102] to describe simultaneously the hardening phenomena in the simple extension of LDPE Melt I with following softening. The first exponential term in (B1) describes the hardening phenomena in polymer melts due to orientation of macromolecules. The simple molecular model and explanations are given in paper [23] and book [32]. The second term in formula (B 1) reflects the softening phenomena in the flow of high oriented melts as described by the Eyring formula in the activated rate processes. As shown in paper [102], the softening is attributed to the thermo-mechanical degradation with a decrease in the molecular mass during extension flow of an LDPE melt. This was revealed by intrinsic viscosity measurements performed on the specimens left over from the experiments. This points to irreversible effects which are outside the scope of a purely rheological description. Interestingly, Wagner had to make a "structural irreversibility"

570

assumption [108] to properly describe recovery following extension flow for LDPE Melt I. To describe simultaneously the flow data for LDPE Melt I, the parameters 13 and v in (B 1) were chosen after a fitting procedure as" [3 - 0.15,

v - 0.03.

(B2)

With the parameters shown in (B2), formula (B1) initially show only hardening and only then, after decay of the first term, it describes "irreversible softening". It should be also noted that in simple shearing, with typical data for shear rates available in rotating instruments, the value b in formula (B 1) is all the time near unity. This fact preserves the description of simple shear data demonstrated in paper [40]. Nevertheless, we manifest in the following also some important comparisons between our calculations and data for simple shearing too. Figure B 1 indicates a good description of the steady shear viscosity and first normal stress difference over a wide range of shear rates. Descriptions of the transient stress growth during shearing experiments are depicted in Figure B2. The discrepancy for N~ at short times is probably due to the instrumental problems for these measurements [97]. 7

I

7

10

~__...~_v

10 6

I

13..

.

~

u

l

~

l

l

10

Melfl o Temp. 150 C -

~eference

o

10

6 ~..,.

o >:, -~

.,..~

10

5

10 5 ~j_~-~7~,, ,--~7_~~^,-, ~

10

c-,-

121

4

104

O

~ "I3 Q

0

m 5=

103

c~

02

10 _

9

10

3

I

u~ ,-,

t~ ol

f-] O

C~ 10

1

.,.

10-4-

I

I

I

I

1

I

I

10-.3

10-2

10 -1

10 ~

101

10 2

10 `3

Reduced

Shear

10 10 4

O ~.~.

R a t e , 7G T (s 1)

Figure B 1 .Steady shear viscosity and the first normal stress coefficient of LDPE Melt I at reference shear rate of l s -1 [97] at various temperatures (symbols).

571

a

10 5

n

c~

i

Melt l Reference

i 0

Temp.

150

C

**t:~:~O

~

z,~

~

~

r"

t..

~

-I

9

o3

.-vaT -

1s

10 2

I0-'

10 0

10 1

10 2

Shear Strain, 7

Figure B2. Transient stress growth normal stress coefficient for LDPE Melt I [97,109]. Various symbols correspond to different temperatures. We now consider shear creep and recovery. The coincidence with experimental data is good for both the shear strain and N1, for the creep condition, as shown in Figure B3. Strain recovery and normal stress decay following unloading after the creep experiment in Figure B3 are shown in Figure B4. Here, the predictions are accurate for strain recovery but N1 is underpredicted. Large amplitude oscillatory shear data, at low to moderate frequencies are indicative of dissipation during cyclical nonlinear deformation. Figure B5 shows a good agreement with the data obtained by Giacomin et al. [ 111 ] for the batch-labeled IUPAC X. We also demonstrate the capability of the model to describe the extension experiments. Figure B6 demonstrates an excellent agreement between calculations and data [112] for extension stress growth, with the use of formulae (B1) and (B2) for simple elongation. Figure B7 shows the comparison between our calculations and data [107] for the stress growth coefficient in planar extension.

572

10

I

10

I

Melt l Temp. 150

Reference

10

o C

4 o "~

12

-

10

6

_

105

_

10 4

Pa

d

z

N1

o-0

10

~" U, 121

lc" (/3

10 0

_

10 -~

10 -2

10

u

i

i

i

10 -1

100

101

102

Reduced

10 3

2

10 5

time, f / a T (s)

Figure B3. Shear strain and first normal stress difference in creep under constant shear stress for LDPE Melt I [110]. Various symbols correspond to experiments at various temperatures. 102

10

.

.

.

105

.

I

I

I

I

1

104 -

N,

t

z

> 0

~ r ~~

10

0

_

10 3

,~ -la Q

r o_ El k_ (z)

10

-1

10 -2 10 -2

Reference Temp. 150

I 10 -1

I

I

I

100

101

102

Reduced

o

C

102

101 103

time, f / a T (s)

Figure B4. Shear strain recovery and first normal stress difference after unloading following the experiment in Figure B3 [110]. Various symbols correspond to various temperatures.

573

80

I

60

-

40

_

I

1

I

I

I

1

70 = 5.0

Frequency = 1Hz

O~

-~

2o

e

0

(o)

r~

-20 o~

-40

IUPAC

X

o

150

-60 I

-80

-40

I

-30

I

-20

-10

C

I

I

I

I

0

10

20

50

40

-t

Shear

100

75 5O v

I

7o

I

=

I

Rate

(s

I

)

I

I

I

10.0 OQ

Frequenc

=

25

f/l I1)

~

0

_

(b)

-25

-

&)

rn

IUPAC X

-50

-

o

150

-75 I

-100

-80

-60

I

-40

I

-20

C

_

I

I

I

I

0

20

40

60

80

-1

Shear

R a t e (s

)

Figure B5. Large amplitude oscillatory shearing of LDPE IUPAC X [ 111 ].

574

10

6I

I

t

I

I

0.I

1.0

I

0.01

I 0 13_

1

..---2 .+r 10 v

+#

f

Melt I (I50 ~

10

I

I

I

I

I

10 -I

100

101

102

103

10 4

Time (s) Figure B6. Tensile stress growth coefficient vs. time for LDPE Melt I [112]. Finally, we also attempted to describe with the use of formulae (B1) and (B2), the data [98] for equi-biaxial elongation for LDPE Melt I. The result of comparison between our calculations and data is shown in Figure B8. It is evident that the calculated curve predicts more hardening in the rheological behavior of the Melt I as compared to the data. The possible reason for that was discussed in the Section 5 of the paper's main text.

575

106

I

0.05 s

-1

0.01 s 0 0

-1

0~ E~ 13_

10

5

.kp

+o_

Melt I

10

(125 ~C)

4

t

100

..

I

101

102

103

Time (s)

Figure B7. Planar tensile stress growth coefficient vs. time for LDPE Melt I [107]. 1.0

07

c

0 "~ c"

0

C: 0

=

0

0.5

x ~

l.a. 12I

~

~o3 IUPAC A (150 ~ 0

.1

,_

'

0.1

0.3

0.5

9

0.7

I

1.0

Strain

Figure B8. Equi-biaxial damping function for LDPE Melt I [98].

3.0