- Email: [email protected]
An elastic-viscoplastic large deformation model and its application to particle filled polymer film
COMPUTATIONAL MATERIALS SCIENCE ELSEVIER
Computational Materials Science 3 (1994) 146-158
An elastic-viscoplastic large deformation model and its application to particle filled polymer film C.G.F. Gerlach, F.P.E. D u n n e * Department of Mechanical Engineering, University of Manchester Institute of Science and Technology, P.O. Box 88, Manchester M60 1QD, UK Received 27 March 1994; accepted 2 May 1994
Abstract
Constitutive equations have been established to describe the behaviour of elastic-viscoplastic polymer materials initially exhibiting amorphous structure. Through the introduction of a state variable representation of back stress arising from network alignment, the constitutive equations model the large deformation behaviour of the material. A semi-automatic, computational technique has been developed for the determination of the material parameters arising in the constitutive equations. The technique requires tensile stress-strain data to be made available. The constitutive equations have been determined for PET film, and inplemented within a finite element solver. The behaviour of particle filled PET film is then investigated using the material model. The results of the computations show that very high levels of stress are generated local to the particle, with applied strain, which result from localised network alignment and strain hardening. The stresses are likely to lead to particle/matrix debonding, or to the nuclation of voids, since the levels of stress quickly become unsustainable by the material, even for normal processing levels of strain. The results also show that the effects of stress relaxation are significant, with a 10% drop in stress occuring under constant strain, over a period equivalent to one quarter of the loading period.
1. Introduction
There are many industries for which it is important to establish accurate modelling techniques for materials processing, since they are essential in the optimisation of the manufacturing operation. This entails maximising efficiencies, minimising cost and delivery time, and yet maintaining product quality. Polymer film production
* Corresponding author.
is an important class of process because of the very large demand for transparent and opaque film for use in packaging, photography and communications applications. The production of polymer film often involves uni-axial and bi-axial stretching at temperatures close to the glass transition, involving true strains of order 2.0. In the case of P E T and P M M A polymer materials, for example, strains in excess of 0.5 in uni-axial drawing near the glass transition t e m p e r a t u r e lead to the development of network alignment, significant strain hardening, and
0927-0256/94/$07.00 © 1994 Elsevier Science B.V. All rights reserved
SSDI 0 9 2 7 - 0 2 5 6 ( 9 4 ) 0 0 0 5 1 -D
C G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
the associated development of anisotropic behaviour. In the case of rigid particle filled film, the particles act as stress raisers leading to localised regions of strain in excess of that seen by the bulk material. This leads to local alignment of the network, and strain hardening, and with further drawing, a progressively increasing local effective stress field. Because of the highly nonlinear stress-strain behaviour of the material at large strains, which shows a rapidly increasing tangent modulus, very large stresses are likely to be developed around the particle, ultimately leading to particle/matrix debonding or void nucleation and growth, or both. The subsequent propagation of material damage and failure is then likely to be influenced by the particle volume fraction, particle clustering, the magnitude and state of the local stress field, and in particular the state of stress tri-axiality, and the network alignment induced anisotropy. The effect of the damage nucleation and growth is sometimes found to lead ultimately to the rupture of polymer film during processing, incurring significant increases in cost and time. It is clearly of interest to manufacturers to establish the important parameters in determining the damage and failure process, and hence first of all, to develop accurate, yet simple, material models to enable the prediction of the complex material behaviour, ultimately enabling the optimisation of the process.
1.1. Available modelling techniques Although considerable research has been carried out for the description of the material behaviour of polymers, the modelling of this behaviour within large-strain finite element code, desirable for a better insight into the manufacturing process, and its control, is considered to be at a rather primitive stage. The majority of constitutive models implemented in large strain finite element code to date are generally phenomenologicai in character (e.g. Mooney [1]), but in addition, fail to model the material's complex primary and tertiary hardening, and fail to account for very significant rate,
147
and relaxation effects. The accurate modelling of these effects is crucial in enabling a proper insight into the manufacture, and its control, of particle filled polymer film. Physically based constitutive models [2-4] that attempt to characterise the complex behaviour of a class of polymer material have been developed, but have not yet found widespread use in industry. Argon [2] developed a micro-mechanical model for polymer plasticity predicting the temperature, pressure and strain rate-dependent intermolecular resistance during primary hardening. Although accurate predictions for the plastic behaviour of hydrocarbon polymers [2] and aromatic polyimides [5] have been presented, an implementation of the model in finite element code has not been reported. A further constitutive model for the large inelastic deformation of glassy polymers has been developed by Boyce, Parks and Argon [3] which extends Argon's theory of polymer yielding using an additional "back stress" term, which accounts for the change in configurational entropy due to the large-scale chain orientations. The model has been used for the simulation of thermo-mechanically couplcd deformation processes for polymers including hydrostatic extrusion processes. The model has been successfully implemented within the commercial code ABAQUS, but requires the determination of ten material parameters, from a range of experimental tests. The use of the constitutive model has been limited to PMMA to date. The three-dimensional, large strain, material model developed by Buckley [4] for the deformation behaviour of PET over the temperature range 50°C to 100°C is physically based, and accounts for the temperature and strain-rate dependency of the yield stress using the Eyring theory of the flow of glassy polymers. Although the large deformation behaviour is correctly predicted using an enhancement of the configurational entropy function developed by Edwards and Vilgis [6], the low strain behaviour, near yield, is not well represented. The model requires the determination of twelve physically based material parameters from a range of experimental tests.
148
C,G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
1.2. Scope of the present work
2.2. Physics of the deformation process
Large strain constitutive equations are established which fully characterise the primary, secondary, and tertiary hardening behaviour of PET polymer film, and which fully account for strain rate effects and stress relaxation. Semi-automatic computational procedures are developed to enable the constitutive equations to be determined from a material data-base containing uni-axial stress strain data over a range of strain rates. The methods developed are general in nature and application, but are employed to determine the equations for PET film using data available in the literature. The constitutive equations are generalised to three dimensions, and implemented within a finite element solver and, once validated, are used to examine the behaviour of low volume fraction particle filled PET film under uni-axial drawing up to strains of 1.2. The relaxation characteristics of the particle filled polymer film are also examined.
In order to derive a physical understanding of the deformation process, the initial spatial arrangement of the isotropic amorphous state may be illustrated as thread-like, linear-chain molecules forming a tangled mass with no order at the local level. Small-scale deformations are then intuitively represented by the uncoiling of the entangled polymer chains in the stretch direction. The uncoiling mechanism can be described physically using the micro-mechanical model developed by Argon [2]. The randomly orientated chain segments are superimposed under loading by molecular kinks producing a small additional chain extension and orientation. The formation of a molecular kink is subjected to the intramolecular resistance of the chain segment and, to the inter-molecular resistance due to the elastic interaction of the chain segment with neighbouring molecules. By expressing the elastic intermolecular interactions as functions of the macroscopic applied stress, the segment rotation or the primary material hardening phenomenon can be predicted successfully. For large-scale deformations however, not only are chain segments extended and rotated, but in addition, chain motion and slippage occurs, and ultimately the polymer structure changes markedly due to network alignment. A randomly orien-
2. The large strain behaviour of PET
2.1. Phenomenological observations A series of uni-axial, constant width, drawing experiments with PET film (Buckley [4]) has shown that PET exhibits the non-linear material uni-axial, stress-strain characteristics typical of amorphous polymers, as shown schematically in Fig. 1. After a short-range elastic response, the material shows true stress-true strain "post-yield" hardening, which stabilises leading to a constant rate of hardening until large strains are achieved. In the tertiary region, significant true strain hardening then takes place, leading to a very rapid increase in stress. No true strain softening has been observed. The onset of significant crystallisation during hot-drawing of PET has been found to coincide with the stretch-induced hardening of the material. At a bulk temperature of 85°C, the onset of tertiary strain-hardening is consistent over the whole range of strain rates at approximately 0.9 true strain.
Tmo Strtut
Fig. 1. Schematic diagram showing the true stress-true strain behaviour of PET film under uni-axial loading with primary, secondary and tertiary regions: I, II and III.
C G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
149
to describe the rate dependent, large strain behaviour of a class of polymer materials. The constitutive equations describe the network alignment induced strain hardening through a state variable representation, in which the network alignment is modelled by a back stress, the evolution of which is described by an additional constitutive relation.
3. Establishment of viscoelastic viscoplastic constitutive equations
N/ Fig. 2. Schematic diagram showing the orientation alignment of the network due to uni-axial loading for (a) the unloaded amorphous state, and (b) the loaded, aligned, anisotropic state. tated molecular network, representing an amorphous polymer, is shown schematically in Fig. 2(a), which is subjected to uni-axial loading, leading to the network alignment shown schematically in Fig. 2(b). The alignment developed leads to significant anisotropy in the resulting material behaviour, and is often accompanied by the formation of crystallites, giving a semi-crystalline material structure. Expressing the deformation mechanism in energy terms, the change in the configurational entropy predominates over the increase in inter-molecular enthalpy. Here, Haward and Thackray [7] developed a physical description which accounts for the configurational entropy by employing the statistical model of the polymer chain network proposed by Wang et al. [8]. The vehicle of the statistical non-Gaussian network theory enabled Wang and then later Haward to express analytically the free energy associated with a textured microstructure. More recently, Boyce et al. [3] generalised Haward's entropic resistance to three dimensions and employed it to predict the strain-induced change in the material properties in the stretch direction and, hence, to represent the hardening phenomenon observed in the material behaviour of polymers. In the following section, simple viscoelastic, viscoplastic constitutive equations are established
3.1. Maxwell f r a m e w o r k
The proposed non-linear constitutive equations are based on the Maxwell principles established in the classical theory of viscoelasticity, although the theory is modified to enable inelastic strains to be included. The total strain is decomposed into an elastic, E~ and viscous, E" component. The viscous component describes both reversible and irreversible time-dependent deformations in a unified manner, hence enabling both reversible viscoelastic deformations (generated through uncoiling of polymer chains, for example) together with irreversible viscoplastic deformations (generated through network alignment and crystallisation) to be modelled: ~=~+~v.
(1)
This implies that the polymer material may be represented with two constituents exposed simultaneously to the same external loading: one constituent behaves purely elastically whilst the other represents the viscoelastic-viscoplastic character. The elastic response, which accounts for the material behaviour observed at very small strains, and which is small compared with the viscous strain, is modelled using the Hookean relationship, written in rate form as: ~c = 6"/E,
(2)
where E is the instantaneous modulus. The reversible viscous material constituent is based on the dashpot model, in which the rate of change of strain is proportional to the applied stress. However, in order to account for the ob-
150
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
served nonlinearity in the stress strain behaviour at low strains, and the creep behaviour of the material under constant stress, a constitutive equation of the following form is adopted: ~v = k i t , t - m ,
(3)
where k, n and m are material parameters without specific physical meaning. This form of equation is introduced, since it has been employed successfully in modelling a wide range of creep behaviour in metals (Hayhurst et al. [9]), and in particular, enables short term strain hardening typical of that observed in PET, to be modelled. However, the constitutive equation is unable to model network alignment and the resulting hardening, which is essential for the modelling of polymer film manufacture, and this will be addressed in the next section. Combining Eqs. (1), (2), and (3) gives i = ( ( r / E ) + k a " t -m.
(4)
The predictions of Eq. (4) qualitatively match very well the stress-strain behaviour of the polymer material over small strains, and in addition, the equation enables the characteristics of stress relaxation to be modelled. 3.2. Modelling o f network alignment-induced strain hardening
able, referred to as the back stress, R (Agah et al. [10]). The progress of continued strain induced orientation is consequently represented through an increase of residual stresses in the polymer structure. The mechanism leading to the evolution of the back stress, described above, indicates that its rate of evolution depends upon its current level. The existence of network alignment in regions defining the crystallites, leads to the preferential development of further alignment in these regions, and its extension to neighbouring regions, hence the evolution of back stress is given as
k a R. (5) In addition, the development of local mismatch in plastic strains between the regions of network alignment and the bulk material depends on the bulk inelastic deformation. Furthermore, the alignment evidently leads to directionality of the back stress, which must therefore be described by a tensor in the general case. The direction of the back stress is determined by the inelastic deformation, and for the quite general case in which the elastic strains are small when compared with the overall strains of order 2.0, and initially for the uniaxial case only, a further condition on the back stress may be expressed as /~a~.
During the deformation of the polymer, the random segment orientations associated with primary hardening (Argon, [5]) establish an inhomogeneous microstructure, featuring randomly distributed, localised strain hardened areas or crystallites. The individual crystallites are exposed to different plastic strains and, in general, to higher levels of plastic strain than the bulk material. Further plastic deformation takes place in the individual crystallites and in the bulk material at different rates, leading to an increasing mismatch of plastic strains throughout the microstructure. From equilibrium conditions, the difference in the plastic strain of individual crystallites must be compensated through different levels of elastic strain or stress respectively. In the design of constitutive equations, the residual-type stress distribution due to the inhomogeneous deformation is modelled by introducing an internal vari-
(6)
Combining Eqs. (5) and (6) leads to (7)
t~ = f l R t ,
where /3 is a constant of proportionality. Eq. (7) may be re-written R = R o exp(/3e),
(8)
where the material constant R 0 indicates the initial orientation of the microstructure of the amorphous polymer, and the material constant/3 accounts for the strain dependency of the alignment and strain-hardening effect.The back stress, R, influences only the inelastic strain, and is incorporated accordingly into the constitutive equation given in Eq. (4) as follows: 6"
= -~ + k ( o ' -
/~ =/3R~.
R ) " t -m,
(9)
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
The equation set in Eq. (9) defines the uni-axial material model for the polymer. In order to gencralise the model to make it applicable over a range of materials, it is necessary to establish techniques by which the six material parameters E, k, n, m,/3, and Ro may be determined from a suitable material data base.
The peak stress rate that occurs is assumed to coincide with the point of highest stress during the tertiary hardening shown in Fig. 1. In addition, the peak stress rate, d t r / d t , in Eq. (10) is given by the minimum of (t7 - R) during tertiary hardening, which occurs when (r = R. If for these circumstances tr is given the label ~, and the strain E, the label g, Eq. (8) may be written In ~ , = l n R o + fl~i
4. Determination of material parameters Three (E,/3, and R 0) of the six material parameters are linked directly to isolated physical phenomena, and hence can be identified accordingly. The remaining constants k, m, n are empirically defined and hence must be identified as such so that the model predictions represent the observed experimental data with a minimum deviation. A semi-automatic identification procedure, which identifies the required experiments, is developed and the procedure is applied to material data for P E T reported in the literature.
151
(11)
in which i indicates that stress and strain take different values depending on the strain rate. Hence, the strain dependency of the back stress, given by/3, and the initial back stress, representing the initial level of alignment (expected to be close to zero) R 0 can then be determined from stress-strain data for the material over a range of strain rates. The optimal values of R 0 a n d / 3 may be determined from a least squares fit of Eq. (11) to the experimental data. 4.2. Determination o f empirical material parameters k, m and n
4.1. Determination o f parameters E / 3 and R o
The elastic material constant E at a particular t e m p e r a t u r e is available from quasistatic measurements of the modulus over a range of temperatures. The material parameters relating to the back stress evolution, R 0 and /3, can be obtained from simple uniaxial tensile test data at constant temperature over a range of strain rates. For a uniaxial stretching process at a constant strain rate, a, the constitutive equations given in Eq. (9) may be rewritten as: do" dt
- E a - E k ( t 7 - R ) ' t -'n
(10)
For a given constant strain rate, Eq. (10) indicates that initially, when t = 0, the stress rate is given by E a . As the stress increases with increasing applied stress, the stress rate drops, until the magnitude of the back stress starts to approach that of the stress, after which the second term on the right hand side of Eq. (10) reduces in magnitude. This leads to an increase in the stress rate corresponding to the tertiary strain hardening.
The remaining material parameters, k, m and n, arc not related to any specific physical phenomena, and must therefore be determined indirectly from the experimental data. In the present work, the parameters are determined from simple uniaxial true stress-true strain data at constant temperature over a range of strain rates. In order to identify the optimal material parameters over the considered strain and strain rate range, a numerical optimisation scheme is developed employing a general purpose N A G optimisation algorithm. An objective function is defined in terms of the experimental, o-e(c,/) and the predicted, o'p(e,j), stresses over a range of strains eii, for a range of strain rates. The subscript i refers to the particular strain level, and j to the particular strain rate. The objective function, O(m, n, k), is defined as q
O(m,n,k)
p
= Y'. Y'~ I t r ~ ( e , j ) - c r p ( e # ) l Z ( y k ) )=li=] (12)
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
152
in which p is the number of strain levels considered at a given strain rate, and q is the total number of strain rates considered. The function Z(rk) is discussed later. The objective function given in Eq. (12) is supplied to the NAG algorithm which procedes to determine rn, n, and k in such a way that O(m, n, k) is minimised. The predicted stresses, %(%), given in Eq. (12) are determined by numerically integrating Fxl. (9). This has been carried out using a fourth order Runge-Kutta integration technique with an automatic time stepping algorithm. Convergence of the optimisation to the correct minimum has been found to depend on the quality of the starting values; poor starting values for the material parameters m, n and k, have been found to lead to convergence to local minima, with inappropriate values for the material paramaters (e.g. m or n < 0), which sometimes in addition, leads to stability problems in the numerical integration. In order to ensure convergence to the global minimum, a further term, Z(Yk), has been introduced into Eq. (12) which provides a contribution to the objective function for unreasonable values of the material parameters. Z(rk) takes the form
where ¥rk<0,
H ( r k ) = 0,
W/k > 0
(14)
in which % is a strain rate dependant weighting factor, and rk (k = 1, 3) are the material parameters m, n and k. Hence, for m, n, k > 0, Z(rk) reduces to the weighting factor, and provides a contribution to O(m, n, k) which simply ensures a balanced optimisation over the whole strain rate range, where necessary. However, for m, n, k < 0, a contribution proportional to m, n and k is provided, in addition. This technique has been found to ensure that the optimised material parameters obtained are physically reasonable. The optimisation procedure is carried out over the entire range of strains, including those for which the strain hardening is significant, and the equations in Eq. (9) are therefore integrated using the values for the material parameters in the back stress evolution equation, /3 and R 0, obtained in the manner described above. In order to provide reasonable starting values for the overall optimisation scheme, initially the technique described above is employed separately for each strain rate, over that range of strains for which strain hardening effects are small. By calculating the average material parameters over all
3
z ( r ~ ) = % + 52 Irk I H(rk),
H(rk)=l,
(13)
k=l
30
True Stress (MPa)
experimental zs
4-4-
predicted (tmi-axial equations, and multi-axial equations implemented
z0
in ABAQUS)
•
0
*
*
I
i
i
I
i
I
i
0.2
0.4
0.6
0.g
1
lJ
1.4
1.6
True Strain Fig. 3. True stress vs. true strain curves for PET obtained at a bulk temperature of 85°C for the strain rates: 0.1 s - 1, 1.0 s - l and 9.0 s - l (a) experimentally (Buckley 1991), (b) predicted by the uni-axial constitutive equations and (c) predicted by the multi-axial constitutive equations implemented in ABAQUS.
C.G.F. Gerlach, F.P.E. Dunne / Computatwnal Materials Science 3 (1994) 146-158
the strain rates, starting values for the material parameters may be determined.
153
All other terms have the meanings described previously. The multiaxial generalisation of the back stress evolution equation is given as
4.3. Application of material model of P E T at 85°C k = 3/3J(R)D Thc material model described above has been used to model the behaviour of P E T film at 85°C over the nominal strain rate range 0.1 s -1 to 9 s-1. The material parameter determination techniques described above have been applied to the uniaxial stress-strain data generated by Buckley [4], and the material parameters determined as: E = 460 MPa, R 0 = 2.92 × 10 -5 MPa,/3 = 10.748, k = 0.064, m = 0.671 and n = 1.222. The comparisons of the stress-strain curves produced by the model, and those determined experimentally by Buckley, are shown in Fig. 3. The computed curves show the ability of the model to correctly predict the strain rate effect, in addition to providing a good representation of the primary hardening, the stabilised secondary hardening, and in particular, the tertiary hardening due to network alignment.
5. Multiaxial generalisation of material model
The rate of deformation tensor, D, is decomposed in terms of thc elastic and viscous constituents as v
(15)
The elastic constituent is given by: V
D e = CF,
(16)
wherc C is the compliance tensor, and /~ is the Jaumann rate of the Kirchhoff stress. The viscoelastic viscoplastic constituent D v is obtained by consideration of a viscous potential 1/' from which the viscous rate of deformation may be determined as D ~= ~ k t - " J ( F - R ) " - I ( F
in which the rate of evolution of the back stress is related to the current level of back stress, for the reasons given earlier, and the directionality of the back stress is defined by the deformation tensor. In this way, the evolution of the back stress is highly strain path dependent, reflecting the anisotropic nature of network alignment. In addition, the different material states generated from non-proportional loadings may also be modelled. 5.2. Implementation within A B A Q U S The material model described above has been implemented within the general purpose nonlinear finite element code ABAQUS, using the user defined subroutine UMAT. In order to ensure stability of the numerical integration of the rate constitutive equations, it was necessary to control the time step size, since an explicit integration scheme was adopted. 5.3. Validation of model
5.1. Constitutice equations
D:D"+D
(18)
' -R')
(17)
in which the (') is used to indicate a deviatoric quantity and JC F - R ) = [ ~( F' - R') : ( F, - R,)] 1/2
The implementation of the material model in A B A Q U S has been subsequently validated against the uniaxial experimental data for PET, and the model predictions made using the uniaxial form of the constitutive equations. The comparison was made by imposing constant, uniaxial •nominal strain rate loading onto a four element, axisymmetric finite element mesh, up to 200% nominal strain. The predicted true stress-true strain curves are superimposed onto those shown in Fig. 3 for the experimental data, and for the uniaxial prediction. G o o d comparisons are achieved.
6. Modelling the behaviour of particle filled polymer film
The behaviour of particle filled polymer film under conditions of uniaxial drawing is of interest
154
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
for several reasons. Firstly, the presence of the particles tends to lead to the development of asperities on the film surface, since the particle diameter is often of the same order as the film thickness. It is important to be able to quantify the asperity shapes because of their influence on the manufacturing process. Secondly, the particles tend to act as sites for the development of material damage through nucleation and growth of voids, or through particle/matrix decohesion. It is important to be able to quantify the stress fields established around the particles in the viscoelastic viscoplastic material during the drawing process, in order to understand the nature of the damage mechanisms discussed above.
272 plane strain, hybrid, reduced integration elements, ensuring the incompressibility requirement is achieved.
6. 2. Results and discussion 6.2.1. Influence of particle on surface geometry On the application of the loading shown schematically in Fig. 4(a), it is expected that an asperity would develop on the top boundary as a result of the presence of the particle. The displacement of the top free boundary as a result of the applied loading is shown in Fig. 5, for a Dis~olocemenf.
d
6.1. Model assumptions, geometry and loading In order to simulate the uniaxial drawing of particle filled PET polymer film, the bulk drawing process has been assumed to be homogeneous and to take place at a uniform, constant, temperature of 85°C. In addition, it is assumed that conditions of plane strain exist such that the strains in the direction normal to the applied loading but within the plane of the film are zero. The volume fraction of glass particles is specified as 0.1% so that it may be assumed that no interaction between particles occurs. In this way, it is possible to consider a single particle only in the model. The model geometry is shown in Fig. 4(a), which through symmetry conditions, shows only one half of the complete model. The glass particle is located close to the top free surface to enable its influence on surface asperity shape to be investigated. The particle is 5 × 10 -6 m in diameter, and its centre is located 7.5 x 10 -6 m below the top free surface. Displacement controlled loading is imposed on the boundary shown with a constant specified displacement rate of 4 0 X 10 - 6 m s - t , as shown schematically in Fig. 4(b). This displacement rate is maintained for 1.6 s, at which point a bulk nominal strain of 1.6 has been achieved. The displacement rate is then put to zero, with the strain level of 1.6 maintained, to give a constant strain hold period. The particle/matrix system is modelled with
Free 1"o~)
~ , ! . m d a r y' ' ' ~
./
~-x
Free 8attain Boundary
PET M o t r l x
Glass Par'title
80
0 0.0
1.0
2.0
3.0
4.0
T'me (s) Fig. 4. (a) Schemantic diagram showing the boundary and loading conditions imposed on the particle filled PET matrix material model. (b) Diagram showing the displacement controlled loading applied as shown in (a).
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
cle are uniform, with the yy stress taking the value 13.4 MPa. The stresses local to the particle can be seen to become very large compared with the bulk, remote stress, with the maximum yy stress approaching 50 MPa. The location of the peak yy stress coincides with that of the peak Von Mises stress, which approaches 40 MPa. The very high levels of stress, compared with the remote stresses, largely result from the localised strain hardening that occurs due to the strains induced by the presence of the particle. The localised high strains lead to alignment of the network, and a corresponding rapid increase in the back stress. Because of the nature of the PET stress-strain characteristics, any further increase in the applied strain leads to rapid increases in the stress levels local to the particle, which very quickly become unsustainable by the material. For the case of a uniform distribution of particles in PET film, it is highly likely, therefore, that within the localities of the highly aligned network, corresponding to the regions of very high levels of stress local to the particles, particle/ matrix debonding or void nucleation occurs as a mechanism for the redistribution of stress. The redistribution of stress may be sufficient to cause stabilisation of the void growth process at the particles, particularly for very low volume fractions of reinforcement, giving rise to a uniform
nominal strain level of 1.6. It can be seen that an asperity of height 2.5 #m is predicted to develop; this height corresponds to half the particle diameter. In addition, Fig. 5 shows that at approximately 15 tzm along the top free boundary, the displacement of the surface remains constant, indicating that the influence of the particle on the surface geometry extends to approximately three diameters, for this particular particle location. Thc results obtained using the viscoelastic viscoplastic model developed in the present paper, are compared in Fig. 5, with those obtained using a well established Mooney-Rivlin rubber elasticity model. The Mooney-Rivlin model employed is described elsewhere [11], but it is clear from Fig. 5 that the predicted surface characteristics appear to be largely independent of the material model employed. However, this is the case only for the prediction of displacement; the stress fields obtained around the particle are discussed in the following section.
6.2.2. Stress" .fields The yy (i.e. the loading direction) stresses and the Von Mises stresses are shown in Fig. 6(a) and (b) for the case of 1.6 nominal strain. This is the point just prior to the strain hold period in the loading. The stress fields remote from the parti-
in x ~ o n (lam)
155
2-5
0..t
01
I 5
I 10
I 15
I 2O
2S
Position along top free boundary(tun) Fig. 5. Graph showing the displacement of the top free boundary of the loaded particle filled PET film predicted by the present model and the Mooney-Rivlin model.
156
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
C.G.F. Gerlach, F.P.E. Dunne / Computational Matertals Science 3 (1994) 146-158
157
Fig. 6. (a) Predicted stresses in the yy direction obtained using the viscoelastic-viscoplastic material model at a nominal strain of 1.60. (b) Predicted Von Mises stresses obtained using the viscoelastic-viscoplastic material model at a nominal strain of 1.60. (c) Predicted stresses in the yy direction obtained using the viscoelastic-viscoplastic material model at a nominal strain of 1.60, after a strain hold period of (1.4 s.
array of voids in the film, the locations of which are determined by the particles. However, higher volume fractions of particles lead to increased levels of triaxiallity, and may lead to void growth and coalescence, ultimately leading to film tearing. The influence of the free boundary on the stress fields can be seen in both Figs. 6(a) and 6(b), and results in the asymmetry of the stress fields, and the generation of slightly higher yy and Von Mises stresses on the boundary.
6.2.3. Stress relaxation during strain holds At the instant at which a nominal strain of 1.6 is achieved, the displacement rate is set to zero,
and the strain level held constant to allow stress relaxation to occur. The loading is shown in Fig. 4(b). The strain hold takes place over a 0.4 s period, and the resulting relaxed stresses are shown in Fig. 6(c). The yy stresses shown are those predicted by the model. Comparison of Figs. 6(a) and 6(c) shows the extent of the relaxation, in which a peak yy stress of 50 MPa relaxes to approximately 44 MPa in 0.4 s. The characteristics of the stress fields appear not to change greatly during the relaxation, except for the magnitudes of the stress. However, Fig. 6(c) shows that the presence of the particle tends to cause the stresses to relax non-uniformly in the regions remote to the particle.
158
C. G. F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
7. Conclusions A large deformation elastic-viscoplastic material model has been established for a class of polymer material, which requires the determination of six material parameters from tensile stress-strain data. Semi-automated computational techniques have been developed for the determination of the parameters. In the development of models for process control in manufacturing, for example, the proposed model offers significant advantages over those of Boyce et al. [3] and Buckley [4] which require the identification of many (ten to twelve) material parameters necessitating the use of complex experimental test programmes. While the present model is applicable only to isothermal processing, constitutive equations of this form have been extended to fully anisothermal loading (Dunne and Hayhurst [12], Dunne et al. [13]), and these methods may be employed for the present model in a similar way. The material model has been implemented in the ABAQUS finite element solver, and used to model the behaviour of particle filled PET film. The predictions show the development of high localised stresses near the particles, that result, in the main, from network alignment and the corresponding increase in back stress. The levels of stress quickly become unsustainable by the material, even under normal processing levels of strain, and it is likely that this leads to particle/matrix debonding and void growth at the particles.
Acknowledgements The authors gratefully acknowledge the contributions to this work of Drs. S.H. Ashdown, P.D.A. Mills, N. Zahlan and D.P. Jones of ICI Films.
References [1] M. Mooney, J. Appl. Phy. 11 (1940) 582. [2] A.S. Argon, Philos. Mag. 28 (1973) 839. [3] M.C. Boyce, D.M. Parks and A.S. Argon, Mech. Mater. 7 (1987) 15. [4] C.P. Buckley, Thermo-viscoelasticity of high performance polymer films during manufacture, Final report on SERC/ICI Co-operative Research Grant GR/E57345. [5] A.S. Argon and M.I. Bessonov, Philos. Mag. 35 (1976) 917. [6] S.F. Edwards and Th. Vilgis, Polymer 27 (1986) 483. [7] R.N. Haward and G. Thackray, Proc. R. Soc. Lond. A. 302 (1967) 453. [8] M.C. Wang and E. Guth, J. Chem. Phy. 20 (1952) 1144. [9] D.R. Hayhurst, P.R. Brown and C.J. Morrison, Proc. R. Soc. Lond., A 311 (1984) 131. [10] A. Agah-Tehrani, E.H. Lee, R.L. Mallett and E.T. Onat, J. Mech. Phys. Solids, 35 (1986) 5, 519. [11] C.G.F. Gerlach and F.P.E. Dunne, Modelling the influence of particle reinforcement on surface geometry in PET film, Internal report DMM.94.17, Dept. of Mech. Eng. UMIST, Manchester. [12] F.P.E. Dunne and D.R. Hayhurst, Proc. R. Soc. Lond., A 437 (1992) 545. [13] F.P.E. Dunne, A.M. Othman, F.R. Hall and D.R. Hayhurst, Int. J. Mech. Sci., 32 (1990) 11,945.
Computational Materials Science 3 (1994) 146-158
An elastic-viscoplastic large deformation model and its application to particle filled polymer film C.G.F. Gerlach, F.P.E. D u n n e * Department of Mechanical Engineering, University of Manchester Institute of Science and Technology, P.O. Box 88, Manchester M60 1QD, UK Received 27 March 1994; accepted 2 May 1994
Abstract
Constitutive equations have been established to describe the behaviour of elastic-viscoplastic polymer materials initially exhibiting amorphous structure. Through the introduction of a state variable representation of back stress arising from network alignment, the constitutive equations model the large deformation behaviour of the material. A semi-automatic, computational technique has been developed for the determination of the material parameters arising in the constitutive equations. The technique requires tensile stress-strain data to be made available. The constitutive equations have been determined for PET film, and inplemented within a finite element solver. The behaviour of particle filled PET film is then investigated using the material model. The results of the computations show that very high levels of stress are generated local to the particle, with applied strain, which result from localised network alignment and strain hardening. The stresses are likely to lead to particle/matrix debonding, or to the nuclation of voids, since the levels of stress quickly become unsustainable by the material, even for normal processing levels of strain. The results also show that the effects of stress relaxation are significant, with a 10% drop in stress occuring under constant strain, over a period equivalent to one quarter of the loading period.
1. Introduction
There are many industries for which it is important to establish accurate modelling techniques for materials processing, since they are essential in the optimisation of the manufacturing operation. This entails maximising efficiencies, minimising cost and delivery time, and yet maintaining product quality. Polymer film production
* Corresponding author.
is an important class of process because of the very large demand for transparent and opaque film for use in packaging, photography and communications applications. The production of polymer film often involves uni-axial and bi-axial stretching at temperatures close to the glass transition, involving true strains of order 2.0. In the case of P E T and P M M A polymer materials, for example, strains in excess of 0.5 in uni-axial drawing near the glass transition t e m p e r a t u r e lead to the development of network alignment, significant strain hardening, and
0927-0256/94/$07.00 © 1994 Elsevier Science B.V. All rights reserved
SSDI 0 9 2 7 - 0 2 5 6 ( 9 4 ) 0 0 0 5 1 -D
C G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
the associated development of anisotropic behaviour. In the case of rigid particle filled film, the particles act as stress raisers leading to localised regions of strain in excess of that seen by the bulk material. This leads to local alignment of the network, and strain hardening, and with further drawing, a progressively increasing local effective stress field. Because of the highly nonlinear stress-strain behaviour of the material at large strains, which shows a rapidly increasing tangent modulus, very large stresses are likely to be developed around the particle, ultimately leading to particle/matrix debonding or void nucleation and growth, or both. The subsequent propagation of material damage and failure is then likely to be influenced by the particle volume fraction, particle clustering, the magnitude and state of the local stress field, and in particular the state of stress tri-axiality, and the network alignment induced anisotropy. The effect of the damage nucleation and growth is sometimes found to lead ultimately to the rupture of polymer film during processing, incurring significant increases in cost and time. It is clearly of interest to manufacturers to establish the important parameters in determining the damage and failure process, and hence first of all, to develop accurate, yet simple, material models to enable the prediction of the complex material behaviour, ultimately enabling the optimisation of the process.
1.1. Available modelling techniques Although considerable research has been carried out for the description of the material behaviour of polymers, the modelling of this behaviour within large-strain finite element code, desirable for a better insight into the manufacturing process, and its control, is considered to be at a rather primitive stage. The majority of constitutive models implemented in large strain finite element code to date are generally phenomenologicai in character (e.g. Mooney [1]), but in addition, fail to model the material's complex primary and tertiary hardening, and fail to account for very significant rate,
147
and relaxation effects. The accurate modelling of these effects is crucial in enabling a proper insight into the manufacture, and its control, of particle filled polymer film. Physically based constitutive models [2-4] that attempt to characterise the complex behaviour of a class of polymer material have been developed, but have not yet found widespread use in industry. Argon [2] developed a micro-mechanical model for polymer plasticity predicting the temperature, pressure and strain rate-dependent intermolecular resistance during primary hardening. Although accurate predictions for the plastic behaviour of hydrocarbon polymers [2] and aromatic polyimides [5] have been presented, an implementation of the model in finite element code has not been reported. A further constitutive model for the large inelastic deformation of glassy polymers has been developed by Boyce, Parks and Argon [3] which extends Argon's theory of polymer yielding using an additional "back stress" term, which accounts for the change in configurational entropy due to the large-scale chain orientations. The model has been used for the simulation of thermo-mechanically couplcd deformation processes for polymers including hydrostatic extrusion processes. The model has been successfully implemented within the commercial code ABAQUS, but requires the determination of ten material parameters, from a range of experimental tests. The use of the constitutive model has been limited to PMMA to date. The three-dimensional, large strain, material model developed by Buckley [4] for the deformation behaviour of PET over the temperature range 50°C to 100°C is physically based, and accounts for the temperature and strain-rate dependency of the yield stress using the Eyring theory of the flow of glassy polymers. Although the large deformation behaviour is correctly predicted using an enhancement of the configurational entropy function developed by Edwards and Vilgis [6], the low strain behaviour, near yield, is not well represented. The model requires the determination of twelve physically based material parameters from a range of experimental tests.
148
C,G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
1.2. Scope of the present work
2.2. Physics of the deformation process
Large strain constitutive equations are established which fully characterise the primary, secondary, and tertiary hardening behaviour of PET polymer film, and which fully account for strain rate effects and stress relaxation. Semi-automatic computational procedures are developed to enable the constitutive equations to be determined from a material data-base containing uni-axial stress strain data over a range of strain rates. The methods developed are general in nature and application, but are employed to determine the equations for PET film using data available in the literature. The constitutive equations are generalised to three dimensions, and implemented within a finite element solver and, once validated, are used to examine the behaviour of low volume fraction particle filled PET film under uni-axial drawing up to strains of 1.2. The relaxation characteristics of the particle filled polymer film are also examined.
In order to derive a physical understanding of the deformation process, the initial spatial arrangement of the isotropic amorphous state may be illustrated as thread-like, linear-chain molecules forming a tangled mass with no order at the local level. Small-scale deformations are then intuitively represented by the uncoiling of the entangled polymer chains in the stretch direction. The uncoiling mechanism can be described physically using the micro-mechanical model developed by Argon [2]. The randomly orientated chain segments are superimposed under loading by molecular kinks producing a small additional chain extension and orientation. The formation of a molecular kink is subjected to the intramolecular resistance of the chain segment and, to the inter-molecular resistance due to the elastic interaction of the chain segment with neighbouring molecules. By expressing the elastic intermolecular interactions as functions of the macroscopic applied stress, the segment rotation or the primary material hardening phenomenon can be predicted successfully. For large-scale deformations however, not only are chain segments extended and rotated, but in addition, chain motion and slippage occurs, and ultimately the polymer structure changes markedly due to network alignment. A randomly orien-
2. The large strain behaviour of PET
2.1. Phenomenological observations A series of uni-axial, constant width, drawing experiments with PET film (Buckley [4]) has shown that PET exhibits the non-linear material uni-axial, stress-strain characteristics typical of amorphous polymers, as shown schematically in Fig. 1. After a short-range elastic response, the material shows true stress-true strain "post-yield" hardening, which stabilises leading to a constant rate of hardening until large strains are achieved. In the tertiary region, significant true strain hardening then takes place, leading to a very rapid increase in stress. No true strain softening has been observed. The onset of significant crystallisation during hot-drawing of PET has been found to coincide with the stretch-induced hardening of the material. At a bulk temperature of 85°C, the onset of tertiary strain-hardening is consistent over the whole range of strain rates at approximately 0.9 true strain.
Tmo Strtut
Fig. 1. Schematic diagram showing the true stress-true strain behaviour of PET film under uni-axial loading with primary, secondary and tertiary regions: I, II and III.
C G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
149
to describe the rate dependent, large strain behaviour of a class of polymer materials. The constitutive equations describe the network alignment induced strain hardening through a state variable representation, in which the network alignment is modelled by a back stress, the evolution of which is described by an additional constitutive relation.
3. Establishment of viscoelastic viscoplastic constitutive equations
N/ Fig. 2. Schematic diagram showing the orientation alignment of the network due to uni-axial loading for (a) the unloaded amorphous state, and (b) the loaded, aligned, anisotropic state. tated molecular network, representing an amorphous polymer, is shown schematically in Fig. 2(a), which is subjected to uni-axial loading, leading to the network alignment shown schematically in Fig. 2(b). The alignment developed leads to significant anisotropy in the resulting material behaviour, and is often accompanied by the formation of crystallites, giving a semi-crystalline material structure. Expressing the deformation mechanism in energy terms, the change in the configurational entropy predominates over the increase in inter-molecular enthalpy. Here, Haward and Thackray [7] developed a physical description which accounts for the configurational entropy by employing the statistical model of the polymer chain network proposed by Wang et al. [8]. The vehicle of the statistical non-Gaussian network theory enabled Wang and then later Haward to express analytically the free energy associated with a textured microstructure. More recently, Boyce et al. [3] generalised Haward's entropic resistance to three dimensions and employed it to predict the strain-induced change in the material properties in the stretch direction and, hence, to represent the hardening phenomenon observed in the material behaviour of polymers. In the following section, simple viscoelastic, viscoplastic constitutive equations are established
3.1. Maxwell f r a m e w o r k
The proposed non-linear constitutive equations are based on the Maxwell principles established in the classical theory of viscoelasticity, although the theory is modified to enable inelastic strains to be included. The total strain is decomposed into an elastic, E~ and viscous, E" component. The viscous component describes both reversible and irreversible time-dependent deformations in a unified manner, hence enabling both reversible viscoelastic deformations (generated through uncoiling of polymer chains, for example) together with irreversible viscoplastic deformations (generated through network alignment and crystallisation) to be modelled: ~=~+~v.
(1)
This implies that the polymer material may be represented with two constituents exposed simultaneously to the same external loading: one constituent behaves purely elastically whilst the other represents the viscoelastic-viscoplastic character. The elastic response, which accounts for the material behaviour observed at very small strains, and which is small compared with the viscous strain, is modelled using the Hookean relationship, written in rate form as: ~c = 6"/E,
(2)
where E is the instantaneous modulus. The reversible viscous material constituent is based on the dashpot model, in which the rate of change of strain is proportional to the applied stress. However, in order to account for the ob-
150
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
served nonlinearity in the stress strain behaviour at low strains, and the creep behaviour of the material under constant stress, a constitutive equation of the following form is adopted: ~v = k i t , t - m ,
(3)
where k, n and m are material parameters without specific physical meaning. This form of equation is introduced, since it has been employed successfully in modelling a wide range of creep behaviour in metals (Hayhurst et al. [9]), and in particular, enables short term strain hardening typical of that observed in PET, to be modelled. However, the constitutive equation is unable to model network alignment and the resulting hardening, which is essential for the modelling of polymer film manufacture, and this will be addressed in the next section. Combining Eqs. (1), (2), and (3) gives i = ( ( r / E ) + k a " t -m.
(4)
The predictions of Eq. (4) qualitatively match very well the stress-strain behaviour of the polymer material over small strains, and in addition, the equation enables the characteristics of stress relaxation to be modelled. 3.2. Modelling o f network alignment-induced strain hardening
able, referred to as the back stress, R (Agah et al. [10]). The progress of continued strain induced orientation is consequently represented through an increase of residual stresses in the polymer structure. The mechanism leading to the evolution of the back stress, described above, indicates that its rate of evolution depends upon its current level. The existence of network alignment in regions defining the crystallites, leads to the preferential development of further alignment in these regions, and its extension to neighbouring regions, hence the evolution of back stress is given as
k a R. (5) In addition, the development of local mismatch in plastic strains between the regions of network alignment and the bulk material depends on the bulk inelastic deformation. Furthermore, the alignment evidently leads to directionality of the back stress, which must therefore be described by a tensor in the general case. The direction of the back stress is determined by the inelastic deformation, and for the quite general case in which the elastic strains are small when compared with the overall strains of order 2.0, and initially for the uniaxial case only, a further condition on the back stress may be expressed as /~a~.
During the deformation of the polymer, the random segment orientations associated with primary hardening (Argon, [5]) establish an inhomogeneous microstructure, featuring randomly distributed, localised strain hardened areas or crystallites. The individual crystallites are exposed to different plastic strains and, in general, to higher levels of plastic strain than the bulk material. Further plastic deformation takes place in the individual crystallites and in the bulk material at different rates, leading to an increasing mismatch of plastic strains throughout the microstructure. From equilibrium conditions, the difference in the plastic strain of individual crystallites must be compensated through different levels of elastic strain or stress respectively. In the design of constitutive equations, the residual-type stress distribution due to the inhomogeneous deformation is modelled by introducing an internal vari-
(6)
Combining Eqs. (5) and (6) leads to (7)
t~ = f l R t ,
where /3 is a constant of proportionality. Eq. (7) may be re-written R = R o exp(/3e),
(8)
where the material constant R 0 indicates the initial orientation of the microstructure of the amorphous polymer, and the material constant/3 accounts for the strain dependency of the alignment and strain-hardening effect.The back stress, R, influences only the inelastic strain, and is incorporated accordingly into the constitutive equation given in Eq. (4) as follows: 6"
= -~ + k ( o ' -
/~ =/3R~.
R ) " t -m,
(9)
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
The equation set in Eq. (9) defines the uni-axial material model for the polymer. In order to gencralise the model to make it applicable over a range of materials, it is necessary to establish techniques by which the six material parameters E, k, n, m,/3, and Ro may be determined from a suitable material data base.
The peak stress rate that occurs is assumed to coincide with the point of highest stress during the tertiary hardening shown in Fig. 1. In addition, the peak stress rate, d t r / d t , in Eq. (10) is given by the minimum of (t7 - R) during tertiary hardening, which occurs when (r = R. If for these circumstances tr is given the label ~, and the strain E, the label g, Eq. (8) may be written In ~ , = l n R o + fl~i
4. Determination of material parameters Three (E,/3, and R 0) of the six material parameters are linked directly to isolated physical phenomena, and hence can be identified accordingly. The remaining constants k, m, n are empirically defined and hence must be identified as such so that the model predictions represent the observed experimental data with a minimum deviation. A semi-automatic identification procedure, which identifies the required experiments, is developed and the procedure is applied to material data for P E T reported in the literature.
151
(11)
in which i indicates that stress and strain take different values depending on the strain rate. Hence, the strain dependency of the back stress, given by/3, and the initial back stress, representing the initial level of alignment (expected to be close to zero) R 0 can then be determined from stress-strain data for the material over a range of strain rates. The optimal values of R 0 a n d / 3 may be determined from a least squares fit of Eq. (11) to the experimental data. 4.2. Determination o f empirical material parameters k, m and n
4.1. Determination o f parameters E / 3 and R o
The elastic material constant E at a particular t e m p e r a t u r e is available from quasistatic measurements of the modulus over a range of temperatures. The material parameters relating to the back stress evolution, R 0 and /3, can be obtained from simple uniaxial tensile test data at constant temperature over a range of strain rates. For a uniaxial stretching process at a constant strain rate, a, the constitutive equations given in Eq. (9) may be rewritten as: do" dt
- E a - E k ( t 7 - R ) ' t -'n
(10)
For a given constant strain rate, Eq. (10) indicates that initially, when t = 0, the stress rate is given by E a . As the stress increases with increasing applied stress, the stress rate drops, until the magnitude of the back stress starts to approach that of the stress, after which the second term on the right hand side of Eq. (10) reduces in magnitude. This leads to an increase in the stress rate corresponding to the tertiary strain hardening.
The remaining material parameters, k, m and n, arc not related to any specific physical phenomena, and must therefore be determined indirectly from the experimental data. In the present work, the parameters are determined from simple uniaxial true stress-true strain data at constant temperature over a range of strain rates. In order to identify the optimal material parameters over the considered strain and strain rate range, a numerical optimisation scheme is developed employing a general purpose N A G optimisation algorithm. An objective function is defined in terms of the experimental, o-e(c,/) and the predicted, o'p(e,j), stresses over a range of strains eii, for a range of strain rates. The subscript i refers to the particular strain level, and j to the particular strain rate. The objective function, O(m, n, k), is defined as q
O(m,n,k)
p
= Y'. Y'~ I t r ~ ( e , j ) - c r p ( e # ) l Z ( y k ) )=li=] (12)
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
152
in which p is the number of strain levels considered at a given strain rate, and q is the total number of strain rates considered. The function Z(rk) is discussed later. The objective function given in Eq. (12) is supplied to the NAG algorithm which procedes to determine rn, n, and k in such a way that O(m, n, k) is minimised. The predicted stresses, %(%), given in Eq. (12) are determined by numerically integrating Fxl. (9). This has been carried out using a fourth order Runge-Kutta integration technique with an automatic time stepping algorithm. Convergence of the optimisation to the correct minimum has been found to depend on the quality of the starting values; poor starting values for the material parameters m, n and k, have been found to lead to convergence to local minima, with inappropriate values for the material paramaters (e.g. m or n < 0), which sometimes in addition, leads to stability problems in the numerical integration. In order to ensure convergence to the global minimum, a further term, Z(Yk), has been introduced into Eq. (12) which provides a contribution to the objective function for unreasonable values of the material parameters. Z(rk) takes the form
where ¥rk<0,
H ( r k ) = 0,
W/k > 0
(14)
in which % is a strain rate dependant weighting factor, and rk (k = 1, 3) are the material parameters m, n and k. Hence, for m, n, k > 0, Z(rk) reduces to the weighting factor, and provides a contribution to O(m, n, k) which simply ensures a balanced optimisation over the whole strain rate range, where necessary. However, for m, n, k < 0, a contribution proportional to m, n and k is provided, in addition. This technique has been found to ensure that the optimised material parameters obtained are physically reasonable. The optimisation procedure is carried out over the entire range of strains, including those for which the strain hardening is significant, and the equations in Eq. (9) are therefore integrated using the values for the material parameters in the back stress evolution equation, /3 and R 0, obtained in the manner described above. In order to provide reasonable starting values for the overall optimisation scheme, initially the technique described above is employed separately for each strain rate, over that range of strains for which strain hardening effects are small. By calculating the average material parameters over all
3
z ( r ~ ) = % + 52 Irk I H(rk),
H(rk)=l,
(13)
k=l
30
True Stress (MPa)
experimental zs
4-4-
predicted (tmi-axial equations, and multi-axial equations implemented
z0
in ABAQUS)
•
0
*
*
I
i
i
I
i
I
i
0.2
0.4
0.6
0.g
1
lJ
1.4
1.6
True Strain Fig. 3. True stress vs. true strain curves for PET obtained at a bulk temperature of 85°C for the strain rates: 0.1 s - 1, 1.0 s - l and 9.0 s - l (a) experimentally (Buckley 1991), (b) predicted by the uni-axial constitutive equations and (c) predicted by the multi-axial constitutive equations implemented in ABAQUS.
C.G.F. Gerlach, F.P.E. Dunne / Computatwnal Materials Science 3 (1994) 146-158
the strain rates, starting values for the material parameters may be determined.
153
All other terms have the meanings described previously. The multiaxial generalisation of the back stress evolution equation is given as
4.3. Application of material model of P E T at 85°C k = 3/3J(R)D Thc material model described above has been used to model the behaviour of P E T film at 85°C over the nominal strain rate range 0.1 s -1 to 9 s-1. The material parameter determination techniques described above have been applied to the uniaxial stress-strain data generated by Buckley [4], and the material parameters determined as: E = 460 MPa, R 0 = 2.92 × 10 -5 MPa,/3 = 10.748, k = 0.064, m = 0.671 and n = 1.222. The comparisons of the stress-strain curves produced by the model, and those determined experimentally by Buckley, are shown in Fig. 3. The computed curves show the ability of the model to correctly predict the strain rate effect, in addition to providing a good representation of the primary hardening, the stabilised secondary hardening, and in particular, the tertiary hardening due to network alignment.
5. Multiaxial generalisation of material model
The rate of deformation tensor, D, is decomposed in terms of thc elastic and viscous constituents as v
(15)
The elastic constituent is given by: V
D e = CF,
(16)
wherc C is the compliance tensor, and /~ is the Jaumann rate of the Kirchhoff stress. The viscoelastic viscoplastic constituent D v is obtained by consideration of a viscous potential 1/' from which the viscous rate of deformation may be determined as D ~= ~ k t - " J ( F - R ) " - I ( F
in which the rate of evolution of the back stress is related to the current level of back stress, for the reasons given earlier, and the directionality of the back stress is defined by the deformation tensor. In this way, the evolution of the back stress is highly strain path dependent, reflecting the anisotropic nature of network alignment. In addition, the different material states generated from non-proportional loadings may also be modelled. 5.2. Implementation within A B A Q U S The material model described above has been implemented within the general purpose nonlinear finite element code ABAQUS, using the user defined subroutine UMAT. In order to ensure stability of the numerical integration of the rate constitutive equations, it was necessary to control the time step size, since an explicit integration scheme was adopted. 5.3. Validation of model
5.1. Constitutice equations
D:D"+D
(18)
' -R')
(17)
in which the (') is used to indicate a deviatoric quantity and JC F - R ) = [ ~( F' - R') : ( F, - R,)] 1/2
The implementation of the material model in A B A Q U S has been subsequently validated against the uniaxial experimental data for PET, and the model predictions made using the uniaxial form of the constitutive equations. The comparison was made by imposing constant, uniaxial •nominal strain rate loading onto a four element, axisymmetric finite element mesh, up to 200% nominal strain. The predicted true stress-true strain curves are superimposed onto those shown in Fig. 3 for the experimental data, and for the uniaxial prediction. G o o d comparisons are achieved.
6. Modelling the behaviour of particle filled polymer film
The behaviour of particle filled polymer film under conditions of uniaxial drawing is of interest
154
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
for several reasons. Firstly, the presence of the particles tends to lead to the development of asperities on the film surface, since the particle diameter is often of the same order as the film thickness. It is important to be able to quantify the asperity shapes because of their influence on the manufacturing process. Secondly, the particles tend to act as sites for the development of material damage through nucleation and growth of voids, or through particle/matrix decohesion. It is important to be able to quantify the stress fields established around the particles in the viscoelastic viscoplastic material during the drawing process, in order to understand the nature of the damage mechanisms discussed above.
272 plane strain, hybrid, reduced integration elements, ensuring the incompressibility requirement is achieved.
6. 2. Results and discussion 6.2.1. Influence of particle on surface geometry On the application of the loading shown schematically in Fig. 4(a), it is expected that an asperity would develop on the top boundary as a result of the presence of the particle. The displacement of the top free boundary as a result of the applied loading is shown in Fig. 5, for a Dis~olocemenf.
d
6.1. Model assumptions, geometry and loading In order to simulate the uniaxial drawing of particle filled PET polymer film, the bulk drawing process has been assumed to be homogeneous and to take place at a uniform, constant, temperature of 85°C. In addition, it is assumed that conditions of plane strain exist such that the strains in the direction normal to the applied loading but within the plane of the film are zero. The volume fraction of glass particles is specified as 0.1% so that it may be assumed that no interaction between particles occurs. In this way, it is possible to consider a single particle only in the model. The model geometry is shown in Fig. 4(a), which through symmetry conditions, shows only one half of the complete model. The glass particle is located close to the top free surface to enable its influence on surface asperity shape to be investigated. The particle is 5 × 10 -6 m in diameter, and its centre is located 7.5 x 10 -6 m below the top free surface. Displacement controlled loading is imposed on the boundary shown with a constant specified displacement rate of 4 0 X 10 - 6 m s - t , as shown schematically in Fig. 4(b). This displacement rate is maintained for 1.6 s, at which point a bulk nominal strain of 1.6 has been achieved. The displacement rate is then put to zero, with the strain level of 1.6 maintained, to give a constant strain hold period. The particle/matrix system is modelled with
Free 1"o~)
~ , ! . m d a r y' ' ' ~
./
~-x
Free 8attain Boundary
PET M o t r l x
Glass Par'title
80
0 0.0
1.0
2.0
3.0
4.0
T'me (s) Fig. 4. (a) Schemantic diagram showing the boundary and loading conditions imposed on the particle filled PET matrix material model. (b) Diagram showing the displacement controlled loading applied as shown in (a).
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
cle are uniform, with the yy stress taking the value 13.4 MPa. The stresses local to the particle can be seen to become very large compared with the bulk, remote stress, with the maximum yy stress approaching 50 MPa. The location of the peak yy stress coincides with that of the peak Von Mises stress, which approaches 40 MPa. The very high levels of stress, compared with the remote stresses, largely result from the localised strain hardening that occurs due to the strains induced by the presence of the particle. The localised high strains lead to alignment of the network, and a corresponding rapid increase in the back stress. Because of the nature of the PET stress-strain characteristics, any further increase in the applied strain leads to rapid increases in the stress levels local to the particle, which very quickly become unsustainable by the material. For the case of a uniform distribution of particles in PET film, it is highly likely, therefore, that within the localities of the highly aligned network, corresponding to the regions of very high levels of stress local to the particles, particle/ matrix debonding or void nucleation occurs as a mechanism for the redistribution of stress. The redistribution of stress may be sufficient to cause stabilisation of the void growth process at the particles, particularly for very low volume fractions of reinforcement, giving rise to a uniform
nominal strain level of 1.6. It can be seen that an asperity of height 2.5 #m is predicted to develop; this height corresponds to half the particle diameter. In addition, Fig. 5 shows that at approximately 15 tzm along the top free boundary, the displacement of the surface remains constant, indicating that the influence of the particle on the surface geometry extends to approximately three diameters, for this particular particle location. Thc results obtained using the viscoelastic viscoplastic model developed in the present paper, are compared in Fig. 5, with those obtained using a well established Mooney-Rivlin rubber elasticity model. The Mooney-Rivlin model employed is described elsewhere [11], but it is clear from Fig. 5 that the predicted surface characteristics appear to be largely independent of the material model employed. However, this is the case only for the prediction of displacement; the stress fields obtained around the particle are discussed in the following section.
6.2.2. Stress" .fields The yy (i.e. the loading direction) stresses and the Von Mises stresses are shown in Fig. 6(a) and (b) for the case of 1.6 nominal strain. This is the point just prior to the strain hold period in the loading. The stress fields remote from the parti-
in x ~ o n (lam)
155
2-5
0..t
01
I 5
I 10
I 15
I 2O
2S
Position along top free boundary(tun) Fig. 5. Graph showing the displacement of the top free boundary of the loaded particle filled PET film predicted by the present model and the Mooney-Rivlin model.
156
C.G.F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
C.G.F. Gerlach, F.P.E. Dunne / Computational Matertals Science 3 (1994) 146-158
157
Fig. 6. (a) Predicted stresses in the yy direction obtained using the viscoelastic-viscoplastic material model at a nominal strain of 1.60. (b) Predicted Von Mises stresses obtained using the viscoelastic-viscoplastic material model at a nominal strain of 1.60. (c) Predicted stresses in the yy direction obtained using the viscoelastic-viscoplastic material model at a nominal strain of 1.60, after a strain hold period of (1.4 s.
array of voids in the film, the locations of which are determined by the particles. However, higher volume fractions of particles lead to increased levels of triaxiallity, and may lead to void growth and coalescence, ultimately leading to film tearing. The influence of the free boundary on the stress fields can be seen in both Figs. 6(a) and 6(b), and results in the asymmetry of the stress fields, and the generation of slightly higher yy and Von Mises stresses on the boundary.
6.2.3. Stress relaxation during strain holds At the instant at which a nominal strain of 1.6 is achieved, the displacement rate is set to zero,
and the strain level held constant to allow stress relaxation to occur. The loading is shown in Fig. 4(b). The strain hold takes place over a 0.4 s period, and the resulting relaxed stresses are shown in Fig. 6(c). The yy stresses shown are those predicted by the model. Comparison of Figs. 6(a) and 6(c) shows the extent of the relaxation, in which a peak yy stress of 50 MPa relaxes to approximately 44 MPa in 0.4 s. The characteristics of the stress fields appear not to change greatly during the relaxation, except for the magnitudes of the stress. However, Fig. 6(c) shows that the presence of the particle tends to cause the stresses to relax non-uniformly in the regions remote to the particle.
158
C. G. F. Gerlach, F.P.E. Dunne / Computational Materials Science 3 (1994) 146-158
7. Conclusions A large deformation elastic-viscoplastic material model has been established for a class of polymer material, which requires the determination of six material parameters from tensile stress-strain data. Semi-automated computational techniques have been developed for the determination of the parameters. In the development of models for process control in manufacturing, for example, the proposed model offers significant advantages over those of Boyce et al. [3] and Buckley [4] which require the identification of many (ten to twelve) material parameters necessitating the use of complex experimental test programmes. While the present model is applicable only to isothermal processing, constitutive equations of this form have been extended to fully anisothermal loading (Dunne and Hayhurst [12], Dunne et al. [13]), and these methods may be employed for the present model in a similar way. The material model has been implemented in the ABAQUS finite element solver, and used to model the behaviour of particle filled PET film. The predictions show the development of high localised stresses near the particles, that result, in the main, from network alignment and the corresponding increase in back stress. The levels of stress quickly become unsustainable by the material, even under normal processing levels of strain, and it is likely that this leads to particle/matrix debonding and void growth at the particles.
Acknowledgements The authors gratefully acknowledge the contributions to this work of Drs. S.H. Ashdown, P.D.A. Mills, N. Zahlan and D.P. Jones of ICI Films.
References [1] M. Mooney, J. Appl. Phy. 11 (1940) 582. [2] A.S. Argon, Philos. Mag. 28 (1973) 839. [3] M.C. Boyce, D.M. Parks and A.S. Argon, Mech. Mater. 7 (1987) 15. [4] C.P. Buckley, Thermo-viscoelasticity of high performance polymer films during manufacture, Final report on SERC/ICI Co-operative Research Grant GR/E57345. [5] A.S. Argon and M.I. Bessonov, Philos. Mag. 35 (1976) 917. [6] S.F. Edwards and Th. Vilgis, Polymer 27 (1986) 483. [7] R.N. Haward and G. Thackray, Proc. R. Soc. Lond. A. 302 (1967) 453. [8] M.C. Wang and E. Guth, J. Chem. Phy. 20 (1952) 1144. [9] D.R. Hayhurst, P.R. Brown and C.J. Morrison, Proc. R. Soc. Lond., A 311 (1984) 131. [10] A. Agah-Tehrani, E.H. Lee, R.L. Mallett and E.T. Onat, J. Mech. Phys. Solids, 35 (1986) 5, 519. [11] C.G.F. Gerlach and F.P.E. Dunne, Modelling the influence of particle reinforcement on surface geometry in PET film, Internal report DMM.94.17, Dept. of Mech. Eng. UMIST, Manchester. [12] F.P.E. Dunne and D.R. Hayhurst, Proc. R. Soc. Lond., A 437 (1992) 545. [13] F.P.E. Dunne, A.M. Othman, F.R. Hall and D.R. Hayhurst, Int. J. Mech. Sci., 32 (1990) 11,945.