A model for alignment between microscopic rods and vorticity

A model for alignment between microscopic rods and vorticity

M. Wilkinson and H. R. Kennard Department of Mathematics and Statistics, The Open University, Milton Keynes, MK7 6AA, England.

Numerical simulations show that microscopic rod-like bodies suspended in a turbulent flow tend to align with the vorticity vector, rather than with the dominant eignevector of the strain-rate tensor. This paper investigates an analytically solvable limit of a model for alignment in a random velocity field with isotropic statistics. The vorticity varies very slowly and the isotropic random flow is equivalent to a pure strain with statistics which are axisymmetric about the direction of the vorticity. We analyse the alignment in a weakly fluctuating uniaxial strain field, as a function of the product of the strain relaxation time and the angular velocity about the vorticity axis. We find that when , the rods are predominantly either perpendicular or parallel to the vorticity.

1 Introduction

Microscopic rod-like bodies suspended in a fluid flow rotate in response to the velocity gradient of the flow. This introduces a degree of order in the orientation of a suspension of particles which can influence its optical or rheological properties. The equation of motion for the orientation of microscopic ellipsoidal particles was obtained by Jeffery [1]. The implications of this equation of motion for the orientation have been considered by numerous authors: for example [2] discusses the motion of general axisymmetric particles, [3] considers the role of Brownian motion, [4, 5] discuss the alignment fields in (respectively) regular and chaotic flows, and [6, 7] are recent experimental contributions which give an extensive list of references. There are, however, still aspects which are not thoroughly understood. One surprising observation (based upon direct numerical simulation (DNS) studies of Navier-Stokes turbulence) is that in isotropic fully-developed turbulence, rod-like particles show significant alignment with the vorticity vector, but negligible alignment with the principal strain axis [8]. This was given a qualitative explanation in [8], but it is desirable to have a model for this surprising effect which can be analysed quantitatively.

This paper considers an exactly solvable model for the alignment of rods with vorticity. The formulation of this model was motivated by observations about the velocity gradient field of turbulence. It has been observed that the fluctuations of the vorticity vector decay much more slowly than fluctuations of the rate of strain: [8] shows evidence that the correlation functions of strain and vorticity both show approximately exponential decay, with decay times and respectively, where is the Kolmogorov timescale of the turbulence. Similar results were reported earlier by Girimaji and Pope [9] and Brunk, Koch and Lion [10]. This observation suggests that it may be helpful to consider the limit as , that is the limit where the vorticity is frozen, in order to explain the observed alignment.

We use an Ornstein-Uhlenbeck process to model fluctuations of the velocity gradient, and consider the limit where the vorticity evolves very slowly. This model is solved exactly in the limit where the strain which occurs over the timescale is small. The alignment of the rod direction and the direction of the vorticity vector can be described by computing the probability density function (PDF) of . We find that in these limits the PDF of , denoted by , can be computed exactly. This analytically solvable model has a single dimensionless parameter, , where is the angular velocity of rotation about the vorticity vector. We find that when , the probability density has two sharp peaks, one at (indicating perfect alignment with vorticity), the other at (implying that the rods are perpendicular to the vorticity). In the limit as , the peak at is higher than at , but it is also narrower, with both peaks containing a finite probability. (Throughout this paper, is the expectation value of , and we use to denote its probability density function).

Section 2 discusses the model which will be solved: the equations of motion for a microscopic rod are considered in section 2.1, and the Ornstein-Uhlenbeck model for the velocity gradient of an isotropic random flow is described in section 2.2. Section 3 discusses a transformation of the equation of motion in which the isotropic velocity gradient is replaced by a pure strain field which is axisymmetric about the direction of the vorticity vector, and it discusses the parametrisation of such axisymmetric random strain fields. Section 4 considers the general solution for alignment of rod-like particles in axisymmetric strain fields, before specialising to the solution of the model developed in section 3. Section 5 summarises our conclusions. The analysis in section 4 is closely related to recent work by Vincenzi [11], who analysed the alignment of ellipsoidal particles in an axisymmetric Kraichnan-Batchelor model.

2 Equations of motion

2.1 Non-linear and linear equations of motion for rods

We consider microscopic objects advected in a fluid with velocity field . The objects are assumed to be neutrally buoyant, and smaller than any lengthscale characterising the fluid, but sufficiently large that their Brownian motion need not be considered. The motion of the body is described by the position of its centre, , and the direction of a unit vector aligned with its axis, . The centre of the body is assumed to be advected by the fluid flow: . The motion of the unit vector defining the axis of symmetry is determined by elements of the velocity gradient tensor, evaluated at the centre of the body:


where is the advected particle trajectory. The equation of motion of the director vector of a microscopic rod-like body is [1]


We assume the flow is incompressible, so that . This tensor can be decomposed into a symmetric part , which is termed the strain rate, and an antisymmetric part , which is the vorticity tensor:


If the velocity gradient matrix were constant in time, the equation of motion (2) would imply that the vector would become aligned with the eigenvector corresponding to its largest eigenvalue. However, numerical simulations of equation (2) for velocity fields of fully developed turbulence show a different, and unexpected, phenomenon [8]. It is found that the direction vector has negligible correlation with the dominant strain eigenvector, but that it does have a quite pronounced correlation with the vorticity vector, .

Our analysis of the alignment due to the motion (2) will use an observation due to Szeri [12]: the non-linear equation (2) can be solved by considering a companion linear equation for a vector , which evolves under the action of a monodromy matrix :


where the initial conditions are (the identity matrix) and . The solution to (2) is obtained by normalising the solution of (4):


The advantage of this approach is that it is easier to solve the linear equation (4) than the non-linear equation (2).

2.2 Ornstein-Uhlenbeck model for velocity gradients in isotropic flows

In this section we describe a simple stochastic model for the matrix in isotropic random flows. A version of this model was used by Vincenzi et al [13], and its structure is suggested by the observations in [10]. The model was also considered in [8], which gave a detailed account of its implementation. Here we give a brief summary.

It is known that the elements of and fluctuate randomly about zero, with different timescales and respectively. Their correlation functions are well approximated by exponential functions. This suggests modelling the elements of and by Ornstein-Uhlenbeck processes [14, 15]. The three independent components of the vorticity will be modelled by:


where the are independent white-noise signals, satisfying


This model predicts that the correlation function of is exponential [14, 15]:


The components of the strain-rate matrix are generated by a further six Ornstein-Uhlenbeck processes, with a different correlation time . The off-diagonal elements are generated by a process of the same form as (6), with the diffusion coefficient in (7) replaced by . The diagonal elements of the strain-rate matrix must satisfy , which is the incompressibility condition, . This constraint is satisfied by the solution of the following Ornstein-Uhlenbeck equations


The elements and generated by these processes are statistically independent, apart from the constraint that . The variances of the off-diagonal, diagonal and vorticity elements are respectively denoted , and , and are related to the relaxation times and diffusion rates by , and . The requirement that the statistics of the model are invariant under rotations (so that it describes a velocity gradient with isotropic statistics) gives , so that this model has four parameters: , , and . Note that the diffusion coefficients have dimension , implying that the model has three independent dimensionless parameters. In the following we consider the limit as , so that the vorticity is frozen, with angular velocity . We also assume that , so that the strain fluctuations are small. This leaves one dimensionless parameter, which we will take to be .

3 Transformation to an axisymmetric pure strain model

3.1 The frozen vorticity limit

In this section we consider the alignment of rod-like particles in an isotropic flow, where there is a non-zero vorticity which is slowly varying. The approach is to transform the equation of motion to a reference frame which rotates around the axis of vorticity. In this coordinate system, the strain field oscillates in directions which are perpendicular to the vorticity vector, in addition to having random temporal fluctuations. The effect of these oscillations is to reduce the effective intensity of the random strain field in directions perpendicular to the vorticity vector, so that an isotropic problem with vorticity is transformed to an axisymmetric model with a velocity gradient which is a pure strain. This reduction was also discussed in [8], but is included here for the convenience of the reader.

In order to isolate the effect of the vorticity in the equation of motion for the monodromy matrix, , we introduce another monodromy matrix which evolves under the vorticity alone:


Note that is just a rotation matrix, describing rotation about an axis in the direction of the vorticity vector . The two monodromy matrices may be related by writing


where is an evolution matrix which describes the effect of the shear. An elementary calculation shows that has the equation of motion




is obtained from by applying a time-dependent rotation. Consider the form of the matrix . In the case where the vorticity vector is frozen, and equal to , the matrix is a rotation matrix: . Without loss of generality we can consider the case where the vorticity is aligned with the -axis, with magnitude , where is the rotational angular velocity, so that is a rotation matrix of the form


If the elements of are , the elements of are


Note that all of the off-diagonal components oscillate with angular frequency or . The diagonal component in the direction of the vorticity vector does not oscillate, but the other diagonal elements contain both oscillatory terms and non-oscillatory terms.

3.2 Limit of short correlation time for strain rate

Now consider the case where the strain rate is sufficiently small that the strain which accumulates over its correlation time is very small. In this case the evolution of the matrix (defined by equation (12)) can be described by a diffusive process. Specifically, we consider the evolution of (12) over a time period which is large compared to the correlation time of the strain field , but still sufficiently small that the strain which accumulates over this time interval is small. We write


where the are small and may be assumed to be random matrices, chosen independently at each timestep. We characterise the evolution (12) by computing the statistics of the random strain increments , which are in turn obtained from the random strain using (12) and (15). The advantage of considering the small elements is that they are small random increments which are applied independently at each timestep. This enables their effect to be analysed using a Fokker-Planck equation. First consider the relation between the elements of the matrices and . By integrating (12) and using the definition (16) we obtain


Iterating this expression, taking , and suppressing the initial time in the arguments of we obtain


Using the fact that the correlation time is assumed to satisfy , we can write


The first of term is a random variable with mean zero and size , giving rise to a diffusion term in a Fokker-Planck equation. The second term represents a drift at a velocity which is well-defined in the limit as . The remaining terms may be neglected. In order to formulate the Fokker-Planck equation, we must determine the statistics of the increments .

If , the effect of the oscillatory terms in equation (15) is negligible. Let us consider how to treat the problem when is not small. To simplify the discussion, consider the quantity


where , and where is a random function which satisfies


The spectral intensity of the fluctuations of is defined by


and we shall assume that , so that . The expectation value of is equal to zero. Its variance is


The third steps assumes that , as well as .

Now consider the effect of the random strain model defined by (6)-(9) in the limit where the timescale of the fluctuations of is very small. We assume that the functional form of the spectral intensity of each component is the same, but that their variances are different, so that the spectral intensity of is , implying that the intensity function is normalised so that . We represent the effect of the randomly fluctuating strain field described by (15) by an effective strain field with diffusive fluctuations. Note that is assumed to satisfy , despite being ‘small’. By applying (23), variance of is


Using the same approach, the full set of non-zero covariances of is


Finally, we must consider the mean values of the increments . As an example, consider the evaluation of . From the second term in the right hand side of (19), we have


Only the diagonal elements of have a non-zero contribution to the mean at : we define velocity coefficients as follows


3.3 Uniaxial random strain in three dimensions

In sections 3.1 and 3.2 we showed how an isotropic model with frozen vorticity and rapidly fluctuating strain can be represented by an axisymmetric model where the velocity gradient is a pure strain . In the limit where the strain which occurs over the correlation time is small, the effect of this strain is represented by a product of matrices , where the small increments are independently distributed at each timestep of size . They have diffusive fluctuations, so that . The matrix is traceless, representing the fact that the velocity field is incompressible. The matrix need not, however, satisfy , although it is clear that the leading order term in (19) is traceless. In this section we discuss how to parametrise such axisymmetric strain fields.

We take this axis of rotational symmetry to be ; the general case is obtained from this one by applying rotation matrices. The strain is described by a matrix , which we write in the form


where , , , and are random variables with mean value zero, and diffusive fluctuations: and , , etc.

Applying a rotation about the axis by angle to the random component of gives a transformed matrix, with elements , , , and , given by


where and . The non-random diagonal component is invariant under rotation about . Note that , so that the element is invariant under rotation.

We require that the statistics of the elements are invariant under the rotation angle . It is clear that and must have the same variance, as must and . Without loss of generality, we can consider a model with . We therefore characterise the model by the following statistics, where , , are three constants:


Other covariances, such as , are equal to zero. The requirement that the statistics of the rotated matrix are independent of leads to the equations


Rotational invariance therefore leads to an equation which must the satisfied by and :


so that the model for a uniaxial random strain has four independent parameters, which we can take to be , , and .

For a special choice of these parameters the model is isotropic. Clearly this requires , and , implying . Also, requiring gives . Solving these equations we find that the the covariances of the random terms are fixed in the isotropic case


Another notable limit of the model is the case where the matrix is diagonal: this model is , implying .

4 Alignment in random strain fields

4.1 General solution in a diffusive axisymmetric strain

In section 3 we described the construction of a model for the alignment of microscopic rods with vorticity, in which the velocity gradient is represented as a strain field with diffusive fluctuations, axisymmetric about the direction of the vorticity. First we consider the alignment of rod-like particles under a succession of independent random shears , which satisfy the conditions derived in section 3.3 for the shear statistics to be uniaxial, before discussing the specific model for rod alignment in section 4.2.

Using the approach summarised by equations (4) and (5), the direction vector of a rod-like particle evolves according to the linear equation


where is the infinitesimal strain in time , previously introduced in equation (17), and is the fractional change in length of the vector under the linear evolution equation. Write , where . Because of rotational symmetry about the -axis, we can assume without loss of generality that the component of is equal to zero. We therefore consider the following orthogonal basis of unit vectors


where is the polar angle, and . Writing , we have


By taking the dot product of (34) in turn with , and , we find, respectively to leading order




Let us characterise the evolution of through the evolution of its projection onto the axis, namely


This is a convenient choice because will have a uniform probability density function for an isotropic strain field. Using (36), we find that


We define the drift velocity and diffusion coefficient of by


Using (42) and (37), (4.1) we obtain




Now consider that statistics of the fluctuations of for the uniaxial strain model. For the model defined in section 3.3, we have


where . We can combine these relations with (42) and (43) to determine and :


where . Using the statistics of the elements , , , and , and ordering the resulting expressions as polynomials in , we have:


The steady-state probability density for , namely , satisfies


In the isotropic case, we have and . In this case we find


and the normalised solution is for .

In the general case, we find that is a factor of , and the differential equation (47) is


it us useful to change the variable to . In terms of , the differential equation (49) may be written


Representing the right-hand-side using partial fractions, we obtain


where are the roots of the denominator on the right-hand-side of (50)


and where the coefficients are


The probability density expressed in terms of is then


where is a normalisation constant.

4.2 Solution of rod alignment model

Now we apply the solution obtained in section 4.1 to the model for alignment of microscopic rods, as developed in sections 2 and 3. In section 2.2 we introduced the Ornstein-Uhlenbeck model for a random, isotropic velocity gradient field. The theory in section 3 made two assumptions. In section 3.1 it was assumed that the vorticity varies slowly, and section 3.2 made a further assumption that the strain field is small. Let us consider the implications of these assumptions for the parameters of the model. The assumption that the vorticity varies slowly implies that is large compared to other timescales in the system of equations. The typical strain rate and the correlation time should satisfy . The solution of the Ornstein-Uhlenbeck process implies


so that the criterion for the strain to be small is simply . The angular velocity is related to the magnitude of the vorticity by . The magnitude of the vorticity is estimated by . The rotation rate has a Gaussian distribution, with variance


Because the Ornstein-Uhlenbeck model has an exponential decay of correlations, given by equation (8), the spectral intensity of the strain fluctuations is a Lorentzian function:


In order to apply the results in section 4.1 we must specify the covariances of the fluctuations of the axisymmetric effective strain tensor. If, in accord with the notation of section 3.3, we normalise the variances so that , , , , the non-zero covariances and expectation values of are


We use the assumption that the original random strain field is isotropic, so that the statistics of these elements satisfy (33). Using (4.2) we obtain