## Abstract

We consider mixing in a viscous fluid by the periodic rotation and translation of a stirrer, the Reynolds number being low enough that the Stokes approximation is valid in the unsteady, two-dimensional flow. Portions of the boundary of the container may also move to contribute to the mixing. The shapes of the stirrer and the container are arbitrary. It is shown that the recently developed *embedding method* for eigenfunction expansions in arbitrary domains is well suited to analyse the mixing properties of such mixers. This application depends crucially on the accurate analytical description of the complex, unsteady field. After carefully validating the proposed method against the recent results found in the literature, examples are given of how the method could be used in practice. A special advantage of the suggested method is that it can be extended to handle three-dimensional mixing flows with virtually no change in the procedure shown here.

## 1. Introduction

A branch of technology where the theory of low Reynolds number flows has wide application is that involving mixing processes. Such processes are widely used in the paint, food, pharmaceutical and paper industries. A great deal of money is spent each year on technology and processes that lead to the better mixing of two or more materials, at least one of which is a viscous fluid. It has been estimated (Harnaby *et al*. 1992; Paul *et al*. 2004) that up to 10 billion dollars yr^{−1} were lost by North American industries alone as a consequence of problems in achieving good mixing. There is, therefore, a practical need to develop methods by which one can predict when and how good mixing can be achieved in low Reynolds number flows.

It would appear that turbulent flows would, naturally, have enhanced mixing properties. This is indeed true and in many situations, for example in gaseous combustors, turbulent flows are preferred. While turbulent flows are obviously more suited to achieving high mixing rates, there are many applications where one is forced to use laminar flows for such purposes. When one is dealing with highly viscous pastes, paints or oils, it would be difficult and energetically expensive to try to achieve turbulent flows for mixing. There are also cases where the high shear rates in turbulent flows would be unacceptable. For example, in some fermentors and bioreactors, the high shear rates in turbulent flows would cause unacceptable cell damage and so a laminar flow regime is essential (Cherry & Papoutsakis 1988). In recent times, the analysis of mixing in microscale devices (e.g. Losey *et al*. 2002; Karniadakis *et al*. 2005) has become important, where the flows are often laminar.

It is well known (Ottino 1989) that chaotic advection is the primary route to achieving the best possible mixing efficiency in the laminar flow regime and that an unsteady advecting velocity field is necessary to produce chaotic advection. The latter requirement is strictly necessary only in two-dimensional flows. Mixing calculations involve (i) the accurate determination of the advecting velocity field, (ii) the calculation of the motion of the particles of the material that has to be mixed, and (iii) an assessment of the mixing characteristics of the flow. That the first step is non-trivial even in the low Reynolds number regime, governed as it is by the simpler Stokes equations, may come as a surprise to non-practitioners; it is often thought that with computational fluid dynamics (CFD) and modern computational resources, the computation of laminar flows should pose no problems. However, in an unsteady mixing flow, the changing field, in a quasi-steady approach, has to be computed repeatedly. Second, the particle-tracking effort can also quickly get out of hand, particularly at the later stages of the mixing process as more particles have to be introduced and tracked; and third, since the velocity field in a CFD approach is in general only numerically known at fixed mesh points, a great deal of computational effort has necessarily to be expended on interpolating the field values at particle locations. It had earlier been estimated that even ‘simple’ mixing flows, in two-dimensional lid-driven containers for which experimental data exist, could take 100 years of computational time on a megaflop machine (Franjione & Ottino 1987). Although the computational situation has improved greatly in recent years, it is still not practically possible to use CFD for reliable mixing calculations. For example, in Zalc *et al*. (2001), Alvarez *et al*. (2002) and Zalc *et al*. (2002), approximately 2 million tetrahedra were used in gridding a stirred tank and the calculation of a single mixing flow took approximately 8 hours each on eight desktop workstations working in parallel. This was for an essentially *steady* three-dimensional field in the rotating frame, where it would have been known in advance that the mixing would be poor. It would be very difficult, if not impossible, at present to compute practical mixing in *unsteady* flows, i.e. ones where the mixing is likely to be good. At present, analytical methods are essential for such work.

It was only recently shown that the eigenfunction expansion method (Shankar 2007*b*), in view of its analytical nature and inherent efficiency and accuracy, could be used to calculate mixing processes, ones for which experimental data were available in Ottino (1989). The cases considered involved viscous fluids in rectangular containers where the essentially two-dimensional flows were generated by the unsteady movements of the lid and bottom wall of the container; the time profiles followed a ‘sin^{2}’ protocol. The mixing patterns obtained from the calculation of the evolution of a small dye blob over 14 periods showed good qualitative agreement with the experimental results given in Ottino (1989). As far as we know, this was the first comparison of this type for rectangular containers and this was only possible because of the efficiency of the eigenfunction method, an analytical one.

However, even this successful effort was only for a mixing process generated by a very simple stirring mechanism. What was considered was mixing in a two-dimensional rectangular container with the fluid motion being generated by the unsteady motion of the container walls. In most practical mixers, the stirring is done by stirrers rotating *inside* the container. These stirred tanks are one of the most widely used pieces of processing equipment (e.g. Tanguy *et al*. 1999; Alvarez *et al*. 2002; Ascanio *et al*. 2005). The stirred tank consists of an impeller or stirrer moving and/or rotating in a tank whose walls may also possibly rotate.

In view of the considerable difficulties involved, there are very few calculations of mixing in stirred containers. Finn & Cox (2001) appear to be the first to successfully deal with the mixing problem in a time-varying domain in a stirred container. They considered a rotating/translating circular stirrer in a circular container full of viscous liquid. Using the method of images in a complex variable formulation of the two-dimensional problem, they were able to give an exact, explicit formula for the quasi-steady field in the stirred liquid. This not only permitted them to calculate detailed mixing patterns of a given initial blob in the fluid, but also to calculate the energetics of the motion. Their work clearly illustrates the advantage, if not the essential need, of having an analytic description of the fluid field. More recently, Cox & Finn (2007) have extended the method to handle paddles of essentially arbitrary cross section, provided only that a conformal mapping of the section to a circle exists.

The calculations described above depend crucially on the external boundary being circular, so that a convenient method of images could be used to satisfy the boundary conditions there; moreover, complex variable-based methods are limited to two-dimensional flows. It is here that the embedding method (Shankar 2005, 2007*b*) used in an eigenfunction expansion formulation may help; all boundaries can, in principle, be arbitrary and the method goes over to three-dimensional flows. However, the technique is not directly applicable to these situations where the domain of the fluid being stirred is not simply connected. Here, the extension to multiply connected domains (Shankar 2006) is needed. With this extension, we are now in a position to consider practical mixing configurations where the stirring is done by stirrers moving slowly, but unsteadily, inside a viscous fluid in a container of essentially arbitrary geometry.

In this paper, we will show how stirring flow fields can be calculated and give a number of examples of mixing calculations based on these. It should be pointed out that the method suggested here is not as explicit or accurate as that of Finn & Cox (2001). However, it has the advantage that the container and stirrer can be of arbitrary shape. We will show that when the container section is circular, we are able to reproduce with good accuracy the results in Finn & Cox (2001). Moreover, the method suggested here should be capable of handling three-dimensional mixing flows.

## 2. The mixing analysis of stirred flow fields

### (a) Formulation and the expansion for the fluid field

Figure 1 is a schematic of the class of flow fields that we will be dealing with. _{0} is the closed boundary of a container that has a liquid of density *ρ* and viscosity *μ* in it. The liquid is stirred by the rotary and translational motion of a stirrer and possibly also by the motion of a part of the boundary of the container. The boundary of the stirrer is denoted by _{1}, while is a scale associated with its angular velocity. The motion is assumed to be planar.

Let be a suitable length scale. We take to be the velocity scale. Dimensional quantities may then be related to the dimensionless quantities as follows:(2.1a)(2.1b)(2.1c)(2.1d)where *p* is the pressure and dimensionless quantities are without tildes. We define the Reynolds number of the motion by(2.2)We will assume that although the motion is unsteady, the Reynolds number is low enough that all the inertial terms in the equations of motion can be considered to be inconsequential. In other words, in this low Reynolds number flow, we will treat the motion as being quasi-steady. In which case, the governing equations are the steady Stokes equations,(2.3a)and(2.3b)Since the motion is planar, it is possible to employ a stream function *ψ*(** x**,

*t*) and do away with the continuity equation (2.3

*a*). The velocity field is then given by

*u*

_{x}(

*x*,

*z*,

*t*)=

*ψ*

_{,z}and

*u*

_{z}(

*x*,

*z*,

*t*)=−

*ψ*

_{,x}, while the stream function satisfies the biharmonic equation(2.4)The problem of determining the stirring field then reduces to determining a biharmonic stream function that instantaneously satisfies the non-penetration and no-slip conditions on the boundaries

_{0}and

_{1}, respectively. Note that in this approximation, the time is felt only through the unsteady boundary conditions. With this understanding, in the interests of brevity, we will in the following stop showing the explicit dependence on time.

The stated linear boundary-value problem is ideal for the application of the embedding method first suggested in Shankar (2005) and then extended to multiply connected domains in Shankar (2006). Let be the boundary of a domain that *contains* the domain bounded by _{0}; let be the inner boundary of an unbounded domain that contains the *exterior* of the domain bounded by _{1}. The suggestion in Shankar (2006) is to use the complete set of *interior eigenfunctions* of the boundary , together with the complete set of *exterior eigenfunctions* of the boundary , to solve the boundary-value problem in the doubly connected domain shown in figure 1. In other words, we write(2.5)where *a*_{n} and *b*_{n}, *n*=1, 2, … are scalars, in general complex, that have to be determined; note that the real parts of the sums have to be taken in (2.5). The scalars are determined by truncating the sums to finite ones and choosing them such that they minimize the boundary data errors on the full boundary _{0}∪_{1}. The minimization procedure is discussed in detail in Shankar (2007*b*). In what follows *N*_{1} and *N*_{2} are the number of terms in the first and second sums of (2.5), respectively, while *M* is the total number of points on the full boundary over which the error is minimized.

For the stirring problem sketched in figure 1, a little more may need to be done. Note that, in this case, (i) the outer boundary is made up of two smooth segments, rather than a single one and (ii) there are corners where the field will lose smoothness. Now if it is possible to choose _{0} itself as and we do so, some care needs to be taken. First of all, in general, since satisfy homogeneous conditions on the boundaries, the velocity induced by these at the corners will always vanish. But the velocity induced by will not in general vanish at these points. And thus, from (2.5), it is clear that there will be difficulties at the corners. It is shown in Shankar (2007*b*) that these can be overcome by using suitable combinations of polynomial solutions of the biharmonic equation. Another difficulty is that, in general, a non-vanishing net flow will occur through these boundaries, which cannot be handled by the eigenfunction expansion. For both of the difficulties to be surmounted, we need to append a field *Ψ*(** x**) to those considered earlier, i.e.(2.6)The form of

*Ψ*(

**), which consists of biharmonic polynomials, and the explicit forms for**

*x**ψ*

^{i}(

**) and**

*x**ψ*

^{e}(

**) in a number of useful coordinate systems are given in Shankar (2007**

*x**a*), which is appended as electronic supplementary material to this paper. Although the expansions in Cartesian and plane polar coordinates are available, for example in Shankar (2007

*b*), it seems that the expansions in elliptic coordinates appear for the first time in Shankar (2007

*a*).

### (b) An example of a stirred flow

A number of examples of various mixer configurations with internal stirrers are discussed in Shankar (2007*a*). These include circular containers, non-circular ones with a smooth boundary and ones with a number of corners; the stirrer geometries include circular and elliptical stirrers and ones with three and four rounded corners. In other words, the examples show that the above suggested procedure works for the kinds of configurations that might arise in engineering practice. Here, we will only consider a single example by way of illustration.

Figure 2 shows an elliptical stirrer of thickness/chord ratio of 1/5, with the pivot offset from the centre of the stirrer, rotating in a rectangular container with a trench-like bottom. Note that the thickness of the stirrer is too small for one to be able to efficiently use circular cylindrical eigenfunctions; here, we use the exterior elliptic coordinate eigenfunction expansion and the Cartesian interior expansion, together with the polynomial expansion. Note that in the case considered, the lid is moving to the right while the stirrer is rotating clockwise to energize the flow in the lower part of the container. The figure shows the field at four angular positions *θ*_{0} of the stirrer; their nature is discussed in the electronic supplementary material. A physical understanding of how and why these flow field changes occur in this time-periodic field may be obtained using the methods suggested in Branicki & Moffatt (2006), who show that local corner analysis often leads to an understanding of the global flow evolution.

Altogether the four frames of the figure show that the streamline patterns change considerably as the stirrer rotates; this means that particles instantaneously on a streamline will be on a different streamline an instant later, suggesting that a mixing action will take place. Of course, only a detailed mixing calculation will be able to tell us how good the action is. In such a calculation, one would use the present method to calculate the field at a large number of angular positions and track particles or blobs used as indicators. The advantage of the present procedure is that this can be done very efficiently. It should also be mentioned that one could, with virtually no extra effort, permit the translation of the pivot so that the stirrer could be both rotating and translating for better mixing. We will later evaluate this method of field calculation by comparing it with those found in the literature.

### (c) Mixing in laminar flows

It is well-known, at least for the last 20 years, that a necessary condition for laminar flows to mix well is that they lead to chaotic advection. Unsteadiness is a necessary condition for such advection in two-dimensional flows, while even steady three-dimensional flows can lead to chaotic advection. Thus, in two dimensions, the flow has to be made unsteady either by moving the container boundaries or by having moving elements, i.e. stirrers, inside the flow. Many laminar mixing studies, theoretical, experimental and computational have been published; Funakoshi (2008) is a recent review of some of these where different quantitative mixing measures are also considered; for a nice discussion of the latter, see also ch. 3 by Szalai, Alvarez and Muzzio in Paul *et al*. (2004).

Although it may appear to be easy to visually tell whether two fluids are mixed well, it turns out to be surprisingly difficult to quantitatively estimate the amount of mixing. In addition, chaos happens to be only a necessary condition for good mixing; there are chaotic flows that do not mix well. Some of the tools used in mixing studies are iterated mappings (IM) or Poincaré sections, Liapunov exponents (LE), location of periodic points especially hyperbolic points and their stable and unstable manifolds, etc. An excellent, practical discussion of various mixing measures (dynamical systems based, statistical and physical) is given in Finn *et al*. (2004) who compared 10 such measures as applied to various mixing protocols in a given stirred flow. They not only show that no single protocol simultaneously optimizes all measures, but also that it would be difficult to prescribe reliable general rules for selecting effective protocols. In their study, they show that three of the measures, IM, LE and mixing variance (MV) correlated well over all the 252 protocols studied. This suggests that the computationally less expensive MV may be used instead of the more expensive IM and LE, at least over their domain and protocols. But this may be more generally true. We will, in this paper, only employ this measure to quantitatively assess mixing.

Any mixing calculation has to begin with a computation of the velocity field. Since mixing calculations often involve tracking fluid particles for long times, having an analytic representation is always advantageous, as no interpolation between fixed grid points needs to be done. We have precisely such representations in the present case as the eigenfunction method of §2*a* yields the instantaneous fields in the form of analytic expansions with the need only to store the expansion coefficients. With the computation of the velocity field in hand, one can calculate the trajectory of an individual particle that is being advected by this velocity field. The coordinates of the particle (*x*, *z*) satisfy the pair of ODEs(2.7)

The solution of (2.7) defines a map *Φ* and is usually called the flow; the solution is formally expressed as ** x**(

*t*)=

*Φ*

_{t}(

*x*_{0}) where

*x*_{0}and

**are vectors locating a particle at the initial and arbitrary times. For standard results from dynamical systems theory that are applicable here, we refer the reader to Funakoshi (2008).**

*x*In summary then, our procedure for studying and estimating the effectiveness of mixing in a stirred flow is as follows. The time period *T*=1 is divided into a number, typically 100, of equal subintervals. For each of these instants, the *N* scalars representing the expansion coefficients in (2.5) or (2.6) are calculated and stored. For every tracer particle that needs to be advected at some time, the system (2.7) is solved using a fourth-order Runge–Kutta procedure with the velocity field given ‘exactly’ by the eigenfunction expansion. Once all the particles have been advected for some time periods, the mixing effectiveness can be computed using some suitable measures. The reason that this procedure is practical and effective is that the field is known analytically and the few hundred real scalars in the expansion at each time step *need to be calculated only once*; at time periods after the first they only need to be recalled from memory.

The limitations, posed by computational errors, in the number of time periods over which reliable calculations can be carried out have to be noted. In any actual calculation, errors can arise from errors in the satisfaction of the boundary conditions, function evaluations, the precision of the arithmetic used, the time integration routine, etc.; most of these errors will arise even when using an ‘exact’ formula for the field. As a result, as time progresses the errors will tend to accumulate and any calculation will tend to progressively become inaccurate. A detailed analysis of this is beyond our scope here. We take a practical position here: we try to keep all calculations as accurate as is feasible, make checks by reducing step sizes or increasing *N*_{1} and *N*_{2}, compare with the results in the literature and, if necessary, delete the few particles that will occasionally hit solid boundaries. As far as we can tell, our results compare very favourably with all the relevant results available in the literature and have a satisfactory internal consistency. In fact, we believe that our results show that at least as far as engineering practice is concerned, where mixing in short time is of principal concern, the proposed method should be a useful tool in both design and analysis.

### (d) Validation of the suggested method

Ideally, we would have liked to validate the suggested procedure by comparing our calculated mixing results with those from an experiment in Ottino (1989), as this source had been used to validate the simpler mixing calculations in Shankar (2007*b*). Unfortunately, this book does not deal anywhere with the more complicated situations considered here with internal stirrers. We have therefore validated our procedure by comparing its results with more recent results that have appeared in the literature of mixing calculations involving internal stirring. In all the calculations below, the domain embedding the container was rectangular, while the finite boundary of the domain embedding the exterior of the stirrer was circular.

Galaktionov *et al*. (1999) considered mixing in a rectangular container of dimensions 2*a*×2*b*, with a centrally located, stationary circular stirrer, where mixing is achieved by the discontinuous alternative but periodic motion of the top and bottom walls of the container; the exact protocol and other details are given in their paper. Their field calculation is based on the method of superposition where use is made of the real, Cartesian eigenfunctions of the biharmonic equation (see Shankar 2007*b*); although the expansion coefficients have to be evaluated numerically, the field is given in terms of an analytic representation that permits mixing calculations to be carried out. In their fig. 9*d* they display, after three periods of motion, the deformation of a circular blob originally centred around a candidate hyperbolic periodic point (0.6383*b*, 0), in their notation. Our result for this mixing calculation is shown in figure 3*a*, which appears to qualitatively agree very well with that of Galaktionov *et al*. (1999).

A much studied problem, in view of the existence of an exact solution to the field, is that of mixing in an eccentric annular mixer (EAM). Here, viscous fluid is contained in a circular cylindrical container and is stirred by the rotation of a stationary, eccentrically located stirrer of circular shape. It is the fact that the external and internal boundaries of the doubly connected domain are circular that permits an exact solution to be written down (Jeffery 1922), which then permits mixing calculations to be performed. Our second comparison is with fig. 2 of Funakoshi (2008), which displays the locations, after five periods of motion, of 3600 tracer particles that were initially placed in the lower half of the fluid domain. Here, both the outer and inner boundaries move discontinuously in time in a manner prescribed in §4.2 of the paper. The results of our calculations are shown in figure 3*b*. The agreement is good, and would have been even better if the exact initial distribution of particles had been known.

A far more stringent test of our method is suggested in Finn & Cox (2001) (F–C in the following) who deal with a translating, rotating mixer (TRM) that is a significant improvement on the EAM. Here, the stirrer, in addition to rotating, translates in the fluid, thereby improving the chances of good mixing in the circular container. Using the method of images in a complex analytic framework, F–C derive an exact, analytic formula for the field in the container and use this to accurately calculate mixing for various mixing protocols. Their results provide an excellent bench mark to test other methods used to calculate mixing in complex, time-varying domains; note that in the previous two examples considered, the mixing domain is fixed in time. In order to compare our results with those of F–C, we have programmed their formulae for the velocity field so that the comparisons can be quantitative. In their fig. 9, F–C give streamline plots in the TRM, corresponding to nine sets of values of certain parameters. We have checked our calculations against all of them and find very good agreement. As an example, figure 4*a* compares the results for case 9 of their table 1. The present results are shown as solid lines, while the results from the formula of F–C are shown as dashed lines; on this scale they are indistinguishable. In their fig. 13, F–C compare the results of their calculation for a particular TRM mixing protocol for the evolution of a dye blob with the results of an experiment that they had carried out. For these parameters, the calculated shape of the dye blob after one period is shown in figure 4*b*. The result using the present method is shown as a solid line, while that from the use of F–C's formula is shown as a dashed line; once again the results are indistinguishable on this scale. The further evolution of the blob after three periods is shown in figure 4*c*, which shows the locations of 10 000 points that were initially equidistant on the circumference of the blob. The result using the F–C formula is shown in figure 4*d*; although figure 4*c*,*d* looks very similar, small differences can now be made out. Overall, however, the results compare very favourably.

It is to be noted that, for all the comparisons above, the container was embedded in a rectangular domain, while the exterior of the stirrer was embedded in a domain with circular inner boundary. Thus, the interior eigenfunctions used were of the Papkovich–Fadle type, while the exterior eigenfunctions were the set of separable solutions to the biharmonic equation in cylindrical polar coordinates. It is of some advantage that the same computational procedure can be used for any reasonable shapes of container and stirrer.

## 3. Examples of the analysis of mixing in complex geometries

We have, in the previous sections, shown how the embedding method, in the framework of eigenfunction expansions, provides an efficient tool to analyse mixing in two-dimensional flows driven at least partly by internal stirrers. We will, in this section, consider two examples where we actually use the method in situations that could well mimic those that arise in engineering practice. In the first, the geometry and kinematics are very simple but the results are counter-intuitive and border on amazing; the second is geometrically more complex and meant only to show the generality of the method. Together, they should demonstrate the scope and power of the suggested method.

### (a) The use of a stator to promote mixing and its optimization

In Shankar (2007*b*, pp. 438–440), an unsteady sin^{2} stirring protocol *C*_{2}, in a rectangular container is considered, for which experimental mixing data are available in Ottino (1989); there is good qualitative agreement between the experimental results after 14 periods and the results of the mixing calculations. However, what was of real interest was that this protocol led to a large island of poor mixing near the top, while a slightly different protocol *C*_{1} led to good mixing everywhere. The exact prescription of the protocols is given in Shankar (2007*b*). The question we wish to address here is: is it possible to make a small change in the geometry so that the mixing would be good with protocol *C*_{2} and would it be possible to ‘optimize’ the resulting mixer?

Figure 5*a*,*b* shows the mixing patterns in the rectangular container under the stirring protocol *C*_{2} discussed above after 5 and 10 periods, respectively. What is shown are the positions of 10 000 points initially distributed, equally spaced, on , shown as a dashed vertical line in figure 5*a*. It is clear that, even after 10 periods, there is a large island near the upper lid in which there is very poor mixing. This figure differs a little from that shown in fig. 16.4*b* of Shankar (2007*b*) only because the initial blob there was different and had no intersection with the island, as is the case here; the picture here is entirely consistent with the Poincaré map shown in the book.

What simple measure can we take to improve the mixing in this situation? Recalling the comparison with the results of Galaktionov *et al*. (1999) discussed in §2*d*, a simple idea is apparent: will a *stationary* cylinder positioned in the container improve matters? Figure 5*c*,*d* shows the mixing patterns that are obtained, with the same initial blob, if a small stator of circular section is placed centrally in the container. The result is astonishing: even after just five periods, the mixing now is far better than the bare case after 10 periods. Now that we have a simple practical way to improve the mixing, two questions of engineering importance immediately arise: (i) can we quantify the measure of mixing that has taken place and (ii) can we optimize the improved mixer by say changing the diameter of the stator or its shape or moving its location?

Before we briefly address the two questions raised above, it is worth looking at the time evolution of the mixing in one particular case shown in figure 6, where the cylinder is off-centre, to emphasize that the location can be anywhere. The figure shows the mixing patterns of 10 000 points on the same initial line *x*=0.8 after one, two, three, …, eight periods of motion. The most striking feature is that considerable mixing has taken place by the end of five periods and that, in the following frames, it is very difficult to see the improvement, especially in (*h*) and (*i*). It is clear that one needs to have a computational measure of the mixing to determine how good or bad things are. The figure also shows that, depending on the application in mind, one might well be satisfied with stopping the procedure after a small number of periods, say five or six. Thus, in practical situations, what is often needed is not a knowledge of the asymptotic behaviour of the system, but the behaviour after a small number of periods.

Since we wish to keep the number of parameters small, we consider only five stators of circular section whose diameters are given in table 1. Note that stator 2 was the one considered earlier in §2*d* when comparisons with the results of Galaktionov *et al*. (1999) were made. In order to quantify the mixing, the whole fluid field was gridded with boxes of size 0.01×0.01; the number of tracer particles in each box from some initial blob, the one made up of 10 000 particles shown in figure 5*a*, was counted as a function of the time period and the MV, , was computed. There is a good argument that, since the fluid domains are different in each case, with less mixing to be done by the stators of bigger diameter, the bigger stators ought to be penalized. A corrected MV, MV_{c}=0.6MV/(0.6−0.25*πd*^{2}), where *d* is the diameter of a stator and 0.6 the bare container area, was also computed. Table 2 shows the MVs in each case after 10 periods of motion. It is straight away clear that all the cylinders have a variance an order of magnitude less than the bare container and of very similar magnitude. It is also clear that stator 3 is the least effective. The table also shows data for stator 2 displaced upwards (2^{u}) and displaced to the left (2^{l}); the latter is in fact the most effective. Note also that the correction for area does push up the biggest stator to a less effective position.

Far more instructive than the table is figure 7, which shows how the MV changes with time in each case. There is a great deal that can be learnt from this data. (i) Except for the bare container, MV_{c} falls monotonically with period. (ii) The most rapid fall in MV_{c} occurs for the largest cylinder, even after the penalization for the reduced mixing volume and this occurs monotonically. (iii) Even for large time, decreasing the cylinder diameter first improves the mixing, then worsens it. However, as *d*→0, the MV_{c} will not go to the bare container value; that is a singular torque at the centre will still improve the mixing, a singular perturbation situation. (iv) If, in an engineering situation, one can be satisfied with a mixing level corresponding to an MV_{c} of say 5 or 6, then one could stop the mixing after five or six periods using any cylinder other than 3. This is an important practical conclusion that comes from this analysis. (v) From figure 7*b*, it appears that the off-centred stators are superior for short times and marginally better at long times. One could investigate this further, but we shall not do so.

### (b) Stirring by an elliptical stirrer in a trench-bottomed container

In the previous example, the geometry considered was very simple and improved mixing was achieved statically. We now consider an example that is geometrically more complicated and where internal stirring is required to achieve good mixing. Consider the geometry considered in §2*b*: mixing in a straight-sided container with a trench-like bottom. In this case, it may not be convenient to have the bottom to be in motion and so a stirrer moving in the container to achieve good mixing would seem a natural option. We therefore consider, as shown in figure 2, an elliptical stirrer rotating about a pivot that, in general, is not coincident with its centre. We now permit, however, the pivot also to be in motion, i.e. the stirrer both rotates and translates in the fluid to achieve mixing. A qualitative feel for the types of unsteady streamline patterns that may occur are indicated in figure 2. We use the same notation as in §2*b* and figure 2, but now allow the pivot coordinates (*x*_{p}, *z*_{p}) to be functions of time; consequently, the pivot and stirrer will have a time-dependent velocity (*U*_{p}, *W*_{p}) in addition to the rotation that will have to be included in the boundary conditions.

Since our brief study here is only to act as an illustration of the power of the method advocated, let us limit the number of parameters considered and deal with the single container and stirrer shown in figure 2, and a single stirrer motion protocol, namely F–C's TRM protocol, eqn (76) of their paper, with the F–C parameters being given by: *n*=4, *r*_{1}=0.14, *r*_{2}=0.1033, *a*_{out}=1 and *T*=1. We also specify the kind of mixing problem that we are interested in: given a linear blob on , what can we do to get ‘good’ mixing in the container? It should be noted that the fluid volume in this case is almost twice that in the previous example and that not much can be done at the bottom of the container. Figure 8*a* shows the positions of the stirrer at equal time intervals over a period, indicating a fairly complex stirring motion. Figure 8*b* shows how the 10 000 particles in the blob are distributed in the container after 10 periods of stirrer rotation and translation, *alone*. There is no question that chaotic advection of the particles has taken place. However, as might well be expected, there are very few particles in the neighbourhood of any of the walls, particularly near the upper corners. Now, a fluid dynamicist would expect poor mixing near stationary walls, since the non-penetration and no-slip conditions would inhibit the advection of particles to such neighbourhoods. However, Gouillart *et al*. (2007) suggested from experiments and a simple model that the effects of the walls are actually felt throughout the flow domain slowing down the whole mixing process. They argue that, although the stirrer while moving through a mixed portion of the fluid, does tend to ‘stretch and fold’ the fluid, wide strips of unmixed fluid from the walls are periodically injected between the fine structures; thus, despite exponential stretching in the bulk these strips, from the stationary walls, cause the MV to decay only algebraically. Some evidence in support of this model can be seen in figure 8*b*, where there is a large island in the interior of the mixed portion along with dark streaks that compromise the mixing process. Comparison with the mixing pattern after five periods, not shown here, suggests that the island will not entirely disappear, even if the process is continued.

So what can one do to improve the situation? From our experience with the previous example and from physical intuition, it would seem that the process might be speeded up if relatively fast-moving fluid were injected into the mixing region around the stirrer. In this case, since the bottom cannot be easily moved, and structurally one would want the sides to be stationary, it would appear that one could try to improve matters by having the lid move in its own plane in order to inject fast-moving fluid into the mixing region. We are therefore led to analysing the effect on the mixer, of a lid velocity *q*_{0}, with all other parameters kept fixed. Figure 9 shows the effect of lid motion on the mixing pattern with the same initial blob as in figure 8; all other parameters have been kept the same. Moreover, in this figure, the patterns seen after five periods are also shown so that one can have some idea of how the mixing progresses with time. One can see that, with *q*_{0}=−0.5, the mixing after 10 periods is already superior to that without any lid motion; although the outline of the island can still be made out, its interior is no longer empty. The dark streaks, indicating regions of excessive particles, have either disappeared or broken up and reduced in intensity. It is also clear that, unlike in the first example considered above, there is considerable improvement between *t*=5 and 10. Although one can see the tendency to fill in the corners and the neighbourhood of the lid even when *t*=5, the mixing is still poor when compared with that seen in figure 8*b*; more time is still clearly needed. Although not shown here, the situation appears to be not as good when the lid speed is less, *q*_{0}=−0.25. Does increasing the lid speed further help? It would appear from figure 9*c*,*d* that there is better mixing when *q*_{0}=−1. Now the island is not clearly visible, even when *t*=5, and there is considerable improvement over the next five periods. The upper left corner has filled, the island has disappeared, but some new dark streaks seem to have surfaced. And one wonders whether further increases in lid speed will improve the mixing. Although not shown here, the situation appears worse when *q*_{0}=−2 or −3. So there appears to be an ‘optimum’ *q*_{0} and only a quantitative analysis of the patterns can guide us in this matter.

Figure 10 shows how the MV decays with time for all the lid speeds considered. The following observations are evident from this figure. (i) Beyond *t*=5, the paddle alone achieves very little extra mixing, in the MV sense. (ii) The lid motion improves the mixing for all the cases considered and, except for the *q*_{0}=−0.25 and *q*_{0}=−2 cases, MV decreases monotonically with *t* over the range considered. (iii) For the three ‘best’ lid speeds, MV is less by an order of magnitude than the bare paddle case at *t*=10. (iv) If low MV by *t*=10 alone were the criterion, it appears that a lid speed of −0.5 is the best; but if uniformity of mixing *over the whole volume* was of importance (and there are ways of estimating this), figure 8 suggests that *q*_{0}=−1 may be ‘better’. (v) That the mixing process is much more complicated than discussed above is evident from the behaviour in the range −1≤*q*_{0}≤−3.

We conclude this example by emphasizing that the choices of container shape, paddle shape and motion, and the parameters used were quite arbitrary; any other reasonable container and paddle shapes could be handled as easily. Whereas we have only considered the use of solid stirrers, it may be of interest to note that Stremler & Chen (2007) showed that, under certain circumstances, stirring can be achieved, in unsteady flow, by the use of suitable passive fluid particles acting as ‘ghost rods’ (see also Clifford & Cox (2006) for the generation of ‘topological chaos’ by suitably placed stationary baffles).

## 4. Conclusion

We have tried in this paper to show that the embedding method, used in the framework of eigenfunction expansions, can be a useful tool in the practical design and analysis of low Reynolds number mixers that employ internal stirrers to achieve good mixing. To be able to perform this function, it was essential that the method be analytical in structure so that the field over a period could be represented in terms of a comparatively small number of coefficients that could then be used repeatedly over many periods; this structure also, crucially, permits the field to be evaluated at the required field point without the need for interpolation. We have also shown with many examples that the method is accurate enough, or can be made accurate enough, to handle chaotic advection over many periods, certainly enough for applications where quick mixing over a few periods is essential. Although we have not discussed this here, much greater accuracy than shown here can be achieved, although naturally with greater cost in the use of more eigenfunctions, smaller time stepping and if necessary greater arithmetic precision. In a practical design process, one could start with relatively crude computations while searching for the most promising candidates, which can then be refined using more eigenfunctions and smaller time steps; this would result in savings of cost and time.

Finally, we conclude by pointing out that the method goes over, in principle, to three-dimensional mixing problems. Although we have so far not looked into the details, the success achieved in Shankar (2007*b*) in handling fully three-dimensional problems gives us reason to believe that three-dimensional mixing problems can be handled by the embedding method. This would require generating appropriate three-dimensional interior and exterior eigenfunctions, but otherwise the method would be very similar. The prospects are intriguing.

## Footnotes

- Received October 23, 2008.
- Accepted December 8, 2008.

- © 2009 The Royal Society