## Abstract

The plastics industry today sees huge wastage through product defects caused by unstable flows during the manufacturing process. In addition, many production lines are throughput-limited by a flow speed threshold above which the process becomes unstable. Therefore, it is critically important to understand the mechanisms behind these instabilities.

In order to investigate the flow of a molten plastic, the first step is a model of the liquid itself, a relation between its current stress and its flow history called a constitutive relation. These are derived in many ways and tested on several benchmark flows, but rarely is the stability of the model used as a criterion for selection. The relationship between the constitutive model and the stability properties of even simple flows is not yet well understood. We show that in one case a small change to the model, which does not affect the steady flow behaviour, entirely removes a known instability. In another, a change that makes a qualitative difference to the steady flow makes only tiny changes to the stability.

The long-term vision of this research is to exactly quantify what are the important properties of a constitutive relation as far as stability is concerned. If we could understand that, not only could very simple stability experiments be used to choose the best constitutive models for a particular material, but our ability to predict and avoid wasteful industrial instabilities would also be vastly improved.

## 1. Introduction

The study of fluid mechanics as we know it today has developed over about three centuries, starting with Newton. Over most of that time, attention has been focused on Newtonian fluids; that is, fluids for which the stress is related to the rate of strain by a simple constant viscosity. This model works very well for fluids whose internal structure is very small and memoryless, such as small-molecule liquids and gases.

Over the last century, industrial processes have increasingly been based on molten plastics, whose internal structure is anything but small. A long molecule, which starts out as a random coil when there is no flow, can be stretched out by flow of the fluid around it. The random Brownian forces acting on it effectively produce an entropic force acting to restore its original shape; this causes a stress in the material that depends on the past history of flow. Modelling a material made of such large molecules and creating a *constitutive relation* between stress and strain history, is a key area of current research.

There are many different ways of modelling complex materials like polymeric fluids. Early models were polarized, either completely empirical (e.g. Bingham 1922; Herschel & Bulkley 1926) or based on analytical principles, which must be followed by any material (e.g. Boltzmann 1874; Jaumann 1905; Oldroyd 1950; Rivlin & Ericksen 1955). A thorough review of the development of these models is given by Tanner & Walters (1998). More recently, much effort has been put into deriving model equations from the microscopic physics of the individual molecules. This gives rise to equations, such as the Pom-Pom (Bishko *et al*. 1997) and Rolie-Poly (Likhtman & Graham 2003) models, which are very closely designed to match the architecture of the molecules, but may be very difficult to solve (or even to reduce to closed form). Typically, all these models are verified by comparing the material's measured behaviour with that of the model in standard flow situations: oscillatory shear flow, shear start-up and (much more difficult to attain experimentally) pure extensional flow.

We have now reached the point where an additional benchmark assessment of the constitutive model is needed, a flow stability problem. As we will see in §2, instabilities are an important problem in polymer processing, and we have recently started to see that matching the steady flow behaviour of a fluid with a constitutive equation is not sufficient to match its stability properties. In §4, we see a case (Wilson & Rallison 1999) in which a constitutive change, which might have seemed minor, completely changes the stability of the flow and in §5, we show new results (Wilson & Fielding 2006) for a system where a constitutive change, which makes a qualitative difference to the steady flow properties, makes only tiny changes to the stability.

## 2. Instabilities

Many plastic products are made by some form of extrusion. The melt that is processed may be extruded either into filaments, which are spun into thin fibres or sheets blown into films, or the extrusion may be an intermediate stage of an injection moulding process. In any of these processes, a smooth, steady flow is required to yield a product with smooth edges, uniform structural strength and a good optical finish. If the flow becomes unstable, the results range from small imperfections to whole production runs having to be scrapped.

For example, plastic bins are made by injection moulding and the solid walls are made using a fountain flow, shown in figure 1*a*. Molten plastic is injected into an empty slit region. As it progresses, the leading edge moves through the mould. The material at the walls is stationary and fresh polymer is injected through the centre of the slit. The polymers which reach the front surface will then move up or down and be laid down on the walls of the mould. There is a stagnation point in the middle of the flow front, separating material which will end up on the top wall from that ending up on the bottom wall. This flow can sometimes exhibit the ‘tiger stripe’ instability in which the stagnation point oscillates up and down (Chang 1994; Grillet *et al*. 2002). The name comes from the resultant finish on the product shown in figure 2.

On a bin, of course, a minor imperfection of finish is not too serious; the bin still works, it was never intended to be beautiful and the product is still fit for market. However, other instabilities can have much more serious consequences. In figure 3, we can see two effects of unstable processing flows. In the experiments of figure 3*a*, carried out by Mavridis & Shroff (1994), a three-layer arrangement of polymeric fluids, as in figure 1*b*, was extruded and then cooled to form a solid film. Where the extrusion flow was stable (figure 3*a*(ii)), the optical quality of the film is almost perfect, whereas the product of an unstable flow (figure 3*a*(i)) may well have to be scrapped. In figure 3*b*, taken from Miller & Rothstein (2004), the surface of the extrudate has undergone the *sharkskin* instability; this optically impaired product, too, is wastage.

## 3. Mathematical stability theory

In this section, we summarize the fundamental principles of linear stability theory in order to apply it to our example systems in the following sections.

Suppose our physical system is governed by a set of coupled nonlinear partial differential equations and boundary conditions(3.1)in which *ϕ* is a vector of physical variables, and the functional can include spatial and time derivatives and boundary conditions, and may be highly coupled. Now, suppose that there is one known solution to this set of equations, which we will call the base state *ϕ*_{0}. Typically, the variables in *ϕ*_{0} are time-independent and have limited dependence on spatial variables.

We now add a tiny perturbation to *ϕ*—in a physical system this might be noise from the surrounding environment; in a computer network, the random loss of data packets—and look for information about its evolution. If all possible small perturbations die away and the system returns to state *ϕ*_{0}, then *ϕ*_{0} is stable.

The essence of linear stability theory is an assumption that the underlying system is regular and the response to a very small perturbation can be deduced by looking only at terms up to linear order in the perturbation. Say we add a perturbation , where *ϵ* is a small parameter. Then (by assumption), we can expand about the base state(3.2)where _{0} is a linear operator acting on , which can depend in an arbitrary way on *ϕ*_{0}. As (*ϕ*_{0})=0, if we neglect terms of order *ϵ*^{2} we are left with a linear system(3.3)In this paper, we will be looking at two-dimensional systems (*x*, *y*) that are unbounded in the *x*-direction (the flow direction) and where the base state *ϕ*_{0} depends only on *y*. In that case, decomposition into Fourier modes leads to the system(3.4)The existence of a solution to this system satisfying both governing equations and boundary conditions is a constraint on *ω* and *k* called the *dispersion relation*. We are using temporal stability, in which we impose a real wavenumber *k* and solve for complex *ω*. If Im(*ω*)>0, the perturbation grows exponentially in time and the system is linearly unstable.

## 4. Blurred interfaces in slit coextrusion flow

Our first example system is a coextrusion flow, used industrially in the manufacture of layered films and printed substrates. We will look at a symmetric three-layer flow for convenience, shown schematically in figure 1*b*; this can be viewed as a two-dimensional model of a core-annular coextrusion or simply as a flow in its own right. The point of interest is not the detail of the flow geometry itself, but the effect that varying the constitutive equation can have on the stability of the flow.

The fluid models we will look at are both based on the Oldroyd-B fluid(4.1)(4.2)(4.3)The variables here are: * U*, the fluid velocity;

*, the stress tensor;*

**S***P*, the pressure;

*, the rate of strain tensor, defined as ; and , the polymer extra stress. The parameters of the model (when inertia is neglected, as here) are:*

**E***η*, the solvent viscosity;

*G*, the elastic modulus; and

*τ*, the polymer relaxation time. The derivative defined by equation (4.3) is the

*upper*-

*convected*tensor derivative, which is the appropriate derivative for line elements, the tensor equivalent of the material derivative.

The Oldroyd-B model was originally derived (Oldroyd 1950) by considering properties that any constitutive equation must satisfy. However, almost uniquely in the field, it can also be derived from a microscopic perspective. If each polymer in a dilute solution is assumed to behave as a pair of (hydrodynamically independent) spheres connected by a linear Hookean spring and with Brownian motion acting on the spheres, the macroscopic stress is given by the Oldroyd-B fluid.

### (a) Instability in pure coextrusion

It is well known (Chen 1991*a*,*b*; Hinch *et al*. 1992; Su & Khomami 1992) that there is an interfacial instability in flows of superposed Oldroyd-B fluids. In symmetric three-layer coextrusion of the type we are considering here, there are five dimensionless parameters: the relative interface position *κ* and, in each fluid, the model parameters *c* and *W*(4.4)in which we have used to denote a global average base-flow shear rate, defined as *U*_{base}(centreline)/*y*_{wall}, and assumed that the two fluids have the same solvent viscosity *η*. We will label fluid 1 as the outer fluid (in the regions *y*_{int}<*y*<*y*_{wall} and −*y*_{wall}<*y*<−*y*_{int}) and fluid 2 as the inner (the region −*y*_{int}<*y*<*y*_{int}).

We will look at a single example system throughout this section. The parameters of this are(4.5)

Owing to the symmetry of the system, we expect our perturbation flows to be either sinuous (*y*-velocity an odd function of *y*) or varicose (an even function of *y*), and for the purposes of this example we will only look at sinuous modes. In the long-wave limit *k*→0, our system is unstable to sinuous perturbations with a growth rate of the order *k*^{2}; in the short-wave limit, it is unstable with a growth rate of the order 1. The growth rate for general wavenumbers is shown as the solid curve in figure 4.

### (b) Variations on the interfacial instability

The instability introduced above has been used as a plausible mechanism for many industrial problems, sometimes rightly and sometimes wrongly. The initial motivation for the comparison we are about to make, which was first published as Wilson & Rallison (1999), came from an attempt to understand the *sharkskin* instability in extrusion, shown in figure 3*b*. Renardy (1987) suggested that the fluid very close to the walls would be depleted in polymer and thus act like a different fluid, so that this interfacial instability could cause the observed irregularities in processing. It has since been shown that the sharkskin phenomenon has its origins on a molecular scale, and is closely related to the adhesion of the polymer on the die wall.

We now look at an analogue to the coextrusion problem above, but one where we only have one fluid, whose properties vary steeply around a nominal interface region. In this way, we hope to explore the behaviour of a more physical system where the ‘interface’ is set up by previous flow and so is not strictly sharp. We use two models representing the extremes of paradigm we could use for such a flow.

In our first model, the fluid is assumed to have material properties which depend only on its original location in the channel. That is, each fluid element has material properties (the parameters *η*, *G* and *τ*) that are assumed to depend in some way on the long-time flow history of that element, but which do not change over the time-scales relevant to our stability analysis. This is a material whose properties respond very slowly to their environment; perhaps structure is built up or destroyed by flow, but over a much longer time-scale than we can investigate.

In the second model, we take the other extreme: the properties of a fluid element respond instantaneously to the flow around it. Rather than a dependence on the initial position in the channel, we now have a direct dependence on the shear rate felt by the fluid. This is an extension of the White–Metzner model (White & Metzner 1963), in which the elastic parameter is a function of shear rate; we have simply added a shear-rate-dependent viscosity term to the model.

It is possible to set up these very different models so that their steady behaviour is identical. Here we will look at a case that is very similar to the standard coextrusion problem, but with a slight blurring of the interface.

If we were to take an (unphysically) infinitely steep variation of the material parameters, then for the first model (slow response to the environment) the system is equivalent to the original coextrusion problem of §4*a*. In figure 4, we show the continuous variation of growth rate as the size of the region over which the properties vary is increased.

In figure 5, we compare the results for our two models. We have chosen parameters that mimic a blurred-interface version of the coextrusion problem with parameters as given in equation (4.5), with the interface spread over a sufficiently small region that we still see instability with the slow-reaction model (as in figure 4). We investigated the whole stability spectrum, that is, all values of Im(*ω*) that are valid for a given wavenumber, *k*. As we can see, the stability properties of the two fluids are similar in many ways; in figure 5*a*, the differences between the stable parts of the spectrum are so small as to be irrelevant to any real flow. In figure 5*b*, we see a startling difference between the fluids; for the fluid with instantaneous response to its surroundings, the interfacial instability has disappeared!

This is a beautiful example of a situation in which a property that is difficult to measure in steady experiments—namely, the rate at which a fluid's material properties respond to the flow—is critical in determining the stability of a relatively simple flow.

## 5. Diffusion in micellar shear-banding flow

Surfactants, which are the active agent in many detergents, are bipolar molecules with hydrophobic head-groups and hydrophilic tail-groups. They tend to migrate to interfaces between different phases of a mixed system: for example, in a washing-up liquid solution, they will initially migrate to the air–water interface at the surface, where they lower surface tension and favour the creation of bubbles. When the water is dirty, they can also migrate to the interface between water and grease, encapsulating the grease and keeping it off the clean plates.

In a clean, closed system with no interfaces, a concentrated solution of these bipolar molecules will form *micelles*, clusters of surfactant molecules whose hydrophobic parts are all in the centre (and so away from the water). These micelles can be spheres, or for higher concentrations long chain-like clusters called *wormlike micelles*.

Experimentally, these solutions commonly exhibit spatially heterogeneous, shear-banded states (Britton & Callaghan 1997). The current thinking is that this behaviour, in which regions of two different shear rates are seen simultaneously in flow of a single fluid, arises from a non-monotonicity in the underlying constitutive relation between the shear stress and shear rate for homogeneous flow (Spenley *et al*. 1993; Makhloufi *et al*. 1995; Mair & Callaghan 1996, 1997).

### (a) Constitutive equations for shear banding

The simplest physical constitutive model1 to mimic this dependence is the Johnson–Segalman (JS) model (Johnson & Segalman 1977). The governing equations(5.1)(5.2)(5.3)are very similar to those we have already seen for the Oldroyd-B fluid, in equations (4.1)–(4.3). The only difference is the change in the derivative of the extra stress tensor, which has changed from the upper-convected derivative of (4.3) to a derivative with a slip parameter, *a*. The limit *a*=1 returns us to the Oldroyd-B model. All the variables have the same meaning here as in the previous section; the new tensor * Ω* is the antisymmetric part of the velocity gradient(5.4)

A sample plot of shear stress against shear rate in simple shear flow for this fluid is given in figure 6. The characteristic non-monotonic form of this curve occurs when *a*^{2}<1 and *Gτ*>8*η*. Homogeneous flow with a shear rate on the decreasing part of the curve (e.g. ) is unstable to one-dimensional perturbations with wavevector in the flow-gradient direction (Yerushalmi *et al*. 1970). The fluid thus separates into a structure comprising bands of differing shear rates and , one on each of the stable, upward-sloping parts of the curve. The interface between the bands has its normal in the flow-gradient direction. The shear stress *T* is uniform across the whole flow, as required by a force balance.

The JS model in its original form contains no mechanism for uniquely selecting the shear stress *T*^{*} at which banding occurs. Instead, in a numerical study of shear banding in such a model, the shear stress in the steady banded state depends strongly on the start-up history (Malkus *et al*. 1990, 1991; Spenley *et al*. 1996; Olmsted *et al*. 2000) and can have a stress anywhere in the range *T*_{1}<*T*^{*}<*T*_{2}, as shown in figure 6. This conflicts notably with experiment, which consistently reveals a highly reproducible banding stress.

For this reason, Olmsted and co-workers (Olmsted & Goldbart 1992; Olmsted & Lu 1997; Olmsted *et al*. 2000; Lu *et al*. 2000) have recently introduced a modification to the standard Johnson–Segalman model, introducing a diffusive *non*-*local* term to the evolution equation (5.2) for ** Σ**(5.5)The diffusion constant is a small parameter; nonetheless, it makes a large qualitative difference to the behaviour of the equations, as we shall see in §5

*b*.

### (b) Steady shear-banded flows: no diffusion

We first look at the original Johnson–Segalman model =0. The base-state flow we shall look at contains just two shear bands: for 0<*y*<*y*_{int}, the flow will be simple shear with shear rate ; and for *y*_{int}<*y*<*y*_{wall}, simple shear with shear rate . The corresponding flow variables are(5.6)(5.7)(5.8)(5.9)in which the first normal stress difference *N* is given by(5.10)The shear stress is *T*^{*} and there is another shear stress condition (from the evolution equation for *Σ*_{xy}), which gives a cubic equation for (5.11)which must be satisfied by both and . Of course, where there are two real roots to a cubic equation, there is a third; but as illustrated in figure 6, the stress–strain relation is decreasing near the middle solution, . This means that, as shown by Yerushalmi *et al*. (1970), a band of this shear rate is unstable to one-dimensional perturbations consisting only of variations in shear rate across the flow, and so this band is never seen in practice.

There is no further information to be gleaned from the Johnson–Segalman fluid; the shear stress is not fixed by the model and can take any value over the range where the cubic equation (5.11) has three real roots.

### (c) Steady shear-banded flows: diffusion

The addition of a small diffusive term in equation (5.5) makes a large difference to the behaviour of the fluid in shear-banding flows. Like the addition of small viscosity to the ideal fluid equations, it raises the order of the highest derivative in the equation set, allowing the possibility of forming boundary layers.

We can still find a steady solution in which all the flow variables are functions only of *y*. Assuming that the shear-rate profile is known, the remaining flow variables are highly similar to those given previously in the absence of diffusion,(5.12)but the equations governing *N*(*y*) and the shear stress become(5.14)(5.15)The two solutions and we had in §5*b* for any value of *T*^{*} are valid in this system, but if we are to see a banded state in which both of these exist within different regions of the flow, then the shear rate must vary in some part of the flow. There is a boundary layer scaling of width (*τ*)^{1/2} over which we expect the derivatives of and *N* to be important; if we introduce a new dimensionless cross-channel variable(5.16)near a nominal interface position, the equations governing and *N* in terms of *ξ* are(5.17)(5.18)This is a fourth-order nonlinear ODE system and the criterion that it must match onto the constant-shear rate bands and (whose values depend on *T*^{*}) as *ξ*→±∞, respectively, is a constraint on *T*^{*}. This constraint was shown (Radulescu & Olmsted 2000) to uniquely determine *T*^{*} for a given set of fluid parameters *η*/*G**τ* and *a*.

In the stability calculations that follow, we will assume here that only one region of each shear rate is formed, and there is then only one interface between high- and low-shear-rate bands, regardless of the presence of diffusion. For definiteness, we assume that there is a band of low shear rate in 0<*y*<*y*_{int}, and a high shear rate band in *y*_{int}<*y*<*y*_{wall}, as in figure 1*c*.

### (d) Linear stability

We saw in the previous sections that the behaviour of shear-banded flows in the Johnson–Segalman fluid is qualitatively changed by the addition of a small diffusive term. In this section, we look at the influence of the same diffusive term on the linear stability of such a flow.

Again, we separate the pure Johnson–Segalman calculation from that with the added diffusion term. The linear system to be solved for the pure Johnson–Segalman fluid is essentially straightforward. Each shear band is treated as a separate fluid with its own base state and the interface between the two is assumed to be a material interface; that is, just as in the coextrusion flow of §4, no fluid parcel may pass from one region to the other, which gives a kinematic boundary condition. This flow shows instability over a wide range of parameter values; the sample system we use here has *a*=0.3, *Gτ*=20*η* and *T*^{*}=0.506158*G*. The shear stress *T*^{*} is chosen to match that which would be selected in the presence of diffusion. For these parameters, the flow is unstable to perturbations of all wavenumbers.

The first investigation of the linear stability problem in the presence of diffusion was a numerical study by Fielding (2005). She found an instability over a finite range of wavenumbers (with the range depending on the diffusion coefficient ) whose growth rate in the limit of small diffusion is markedly similar to that of the interfacial instability of a Johnson–Segalman fluid; one set of her numerical results, along with the pure Johnson–Segalman growth rate, is shown in figure 7.

This equivalence between the singular limit →0 and the pure system =0 is surprising, precisely because the limit is a singular one. The coefficient multiplies the highest derivatives in the linear stability system, just as it does in the nonlinear equations governing the base state. Yet, in the base state, the presence of an infinitesimal non-zero diffusion makes a qualitative difference to the response of the system, whereas the stability quantities appear unaffected.

In fact, very recently, we were able to prove this equivalence (Wilson & Fielding 2006) as →0 to leading order in , for wavenumbers *k* for which *k*^{−1}≥(*τ*)^{1/2}. The details are given in full in that work, but the basic idea is laid out below. The linear equations to be solved for the perturbation variables (in lower case) are(5.19)(5.20)(5.21)We split the solution into two parts: an advected term proportional to the *y*-derivative of the base-state and the remainder (denoted here with a tilde). Thus(5.22)(5.23)The coefficient *ζ* is chosen to be independent of *y* but, like the other perturbation variables, it is proportional to exp[i*kx*−i*ωt*].

The resulting equations are(5.24)(5.25)(5.26)Now, the only appearance of the steep derivatives of base-state quantities is on the right-hand side of equations (5.24) and (5.26).

In the base state, there are two large regions of the flow where , *N* and * Σ* are constant. In these regions, most of the right-hand side forcing terms are 0, and if we assume that the term ∇

^{2}

*σ*has only a weak effect, the governing equations are close to those for a Johnson–Segalman fluid without diffusion.

Between these two regions is a narrow layer around *y*_{int} in which the base state quantities vary quickly and the right-hand side driving terms are important. If we make two assumptions: that does not vary quickly within this layer, and none of the variables is singular, we can set *ζ* to eliminate the driving term of (5.26) at the centre of this layer(5.27)and integrate (5.24) across this narrow layer to obtain the conditions(5.28)This means that we can treat our new system (to leading order, as we have neglected the term ) as a layered system of two Johnson–Segalman shear bands with the interfacial boundary conditions being the continuity of(5.29)which are the continuity of total velocity and traction across a boundary located at *y*=*y*_{int}+*ζ* exp[i*kx*−i*ωt*]. The variable *ζ* is defined by equation (5.27), which reduces to the kinematic definition of displacement of a material interface in the limit *k*^{2}→0. Thus, in the low-diffusion limit the leading-order stability system is exactly equivalent to a no-diffusion system with a material interface between the two shear bands.

This is a surprising result simply because the addition of diffusion made such a marked difference to the base-state behaviour of the flow. Two seemingly disparate systems have been shown to have exactly equivalent stability properties.

## 6. Discussion

We have shown two simple examples of the linear stability properties of flows of viscoelastic fluids. In each example, there are two different ways to model the same basic flow. For the coextrusion problem of §4, two fluids whose base-flow properties are identical have completely different stability properties; a fluid whose material properties have a long memory can show instability where a similar fluid with a shorter memory will not. In the shear-banding flow of §5, on the other hand, the addition of stress diffusion, which qualitatively alters the steady-shear behaviour of the fluid, is shown to make only the smallest quantitative change to the stability properties of the flow.

In both cases, the stability properties of the flow could not have been predicted just from the steady behaviour. Stability calculations sample the fluid properties in a different way.

There are several benchmark flows to which any new constitutive equation is subjected: steady and oscillatory shear, pure extension, and so on. I would argue that the stability properties of the fluid in steady flows should be added to this list. We would have an extra set of data available when trying to choose a constitutive model to fit a given fluid, but there are other benefits of just as great an importance. If we were to gain a thorough understanding of the stability properties of all the myriad of constitutive equations in use today, we would be well equipped to combat costly and wasteful industrial instabilities.

An understanding is already developing of how molecular structure relates to constitutive modelling. This allows manufacturers, to a limited extent, to tailor the constituent materials they use to the desired qualities of the finished product. If we spend time studying the relationship between the constitutive equation and stability, we will create a body of knowledge that can be linked to this existing understanding. The result will be a correspondence between molecular structures and stability; information telling us which molecular shapes will work well in different processing flows. In turn, industries can then use *process stability* as one of the criteria for selecting materials, or redesign each process to avoid the unstable flows. Know your enemy: the first step in saving the mass of waste products in the plastics industry is to understand exactly what is causing the processing flows to fail.

## Footnotes

One contribution of 23 to a Triennial Issue ‘Mathematics and physics’.

The electronic supplementary material is available at http://dx.doi.org/10.1098/rsta.2006.1892 or via http://www.journals.royalsoc.ac.uk.

↵By physical in this context, we mean a model that obeys the principle of Material Frame Indifference.

- © 2006 The Royal Society