ISRN Astronomy and Astrophysics

Volume 2012 (2012), Article ID 987275, 10 pages

http://dx.doi.org/10.5402/2012/987275

## On Estimating Fluxes due to Small-Scale Turbulent Convection in a Rotating Star

Institute of Astronomy & Department of Applied Mathematics and Theoretical Physics, University of Cambridge, Madingley Road, Cambridge CB3 0HA, UK

Received 30 November 2011; Accepted 15 January 2012

Academic Editors: R. Avila and S. Bogovalov

Copyright © 2012 D. O. Gough. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

The way in which turbulent fluxes are usually represented in computations of large-scale flow in the convection zones of the sun and other stars is briefly described. A model of an ensemble of eddies that is capable of generalization to circumstances more complicated than the usual essentially spherically symmetrical convection zone is outlined. Generalization usually requires the introduction of new postulates, and, in so doing, also lays bare some of the assumptions, often implicit, in the usual mixing-length formalisms.

#### 1. Introduction

Computations of large-scale flow in stellar convection zones, be it a meridional circulation or simply rotation, require some representation of the fluxes of heat and momentum by turbulent motion on scales too small to be resolved. It is not uncommon in astrophysics to adopt the attitude that since the theory upon which any such representation is based (usually mixing-length theory) cannot be trusted, it is hardly necessary to exercise care in deriving the formulae for the fluxes. Readers can thus be faced with the prospect of learning of the results of complicated numerical computations, yet being unable to know precisely what those results mean.

The main point I wish to make is that the situation can be improved. The so-called theory used for calculating the turbulent fluxes can be refined, and the numerical machinery that has been built to model the large-scale flow can be used experimentally to test the refinements.

First I shall describe the broad principles behind many of the computations of turbulent fluxes in astrophysics. Then I shall outline a procedure whereby the representations of the small-scale motion might be developed in a coherent manner.

#### 2. The Turbulent Fluxes

The first step in all computations of large-scale astrophysical flow, whether it is stated explicitly or not, is to imagine a separation of scales. Motion with characteristic length and time scales greater than some value, usually implicitly prescribed by the numerical technique, is treated by solving the governing equations directly, using some parametrized approximation to describe the motion with characteristic scales smaller than this value: ignoring the small-scale motion completely is a trivial form of such a parametrization. The large-scale fields of velocity and temperature, say, may be considered to derive from the total fields and via the application of an averaging procedure, which I denote by an overbar. Thus I shall refer to fields such as and as mean quantities; the fields and are fluctuations. In this account, I have in mind particularly the convective envelopes of cool stars, in which turbulent motion is driven predominantly by buoyancy.

For simplicity of presentation I refer the fields to rectangular Cartesian coordinates , with vertical, and I shall ignore the curvature of the convection zone. This is a fair approximation if the scales to be parametrized are small enough. Equations of motion may then be written, in an obvious notation, where and are density and pressure, and are specific internal energy and enthalpy (viz., internal energy and enthalpy per unit mass), is the magnitude of the gravitational acceleration , is the radiative heat flux, is time, and is the Kronecker delta. Throughout my discussion I adopt the summation convention for repeated indices except where one of them is enclosed in parentheses. The total energy equation (3) is derived by combining the kinetic-energy equation, obtained by taking the (scalar) product of (2) with to relate the rate of change of the kinetic energy to the rate of working, with the thermal-energy equation which is a statement of the first law of thermodynamics. In (5), is the specific heat at constant volume. These equations must, of course, be supplemented by an equation of state and an equation determining . Viscous terms have been ignored; so too have perturbations to .

To obtain equations for the large-scale flow the averaging procedure is applied to the full equations of motion. Usually some approximations are made in representing the fluctuations, the most common being the neglect of density fluctuations, except where they appear in the enthalpy flux. This constitutes part of the Boussinesq approximation, which has been justified by Spiegel and Veronis [1] for the case in which , when the fluctuations are generated by buoyancy alone, and when the length scales of the fluctuations are all much smaller than the scale heights of density and pressure of the averaged state. Despite the fact that these conditions are believed not to be satisfied in most stars, I too shall nevertheless adopt the practice, obtaining The only fluxes arising from the small-scale turbulent motion that need to be calculated for closing these equations are evidently the Reynolds stress , the turbulent heat flux , where the prime denotes a turbulent fluctuation, and the turbulent kinetic-energy flux .

In the absence of a large-scale mean flow the equations for the fluctuations are given in Boussinesq approximation by [1] where is the specific heat at constant pressure and is the superadiabatic lapse rate; is the inverse dimensionless isothermal compressibility. In what follows I shall assume that the perturbations are optically thick, so, ignoring, for simplicity, fluctuations in the thermal conductivity, where , in which is the first radiation constant, is the speed of light, and is the Rosseland mean opacity. Various procedures have been suggested for treating optically thin perturbations (e.g., [2–5]); I do not discuss them here other than to suggest the procedure recommended by Gough [6], which is consistent with the procedure adopted here for dealing with the other fluctuation equations and which, unlike the other procedures, provides a smooth transition from optically thick to optically thin fluctuations.

#### 3. Rough Estimates of the Turbulent Fluxes

In most cases some kind of diffusion approximation based on mixing-length theory is used to compute the turbulent fluxes. Let us consider first the convective heat flux. Since the Boussinesq approximation has already been used to derive (6), one might as well continue doing so and, in particular, neglect the contribution to the convective heat flux from pressure fluctuations. That flux can then be written (It is perhaps worth pointing out that for a perfect unionizing gas the pressure fluctuation makes precisely no contribution to when the latter is expressed in terms of the temperature fluctuation.) In this and all subsequent equations overbars are omitted from mean quantities, except where there is risk of confusion.

If there were no large-scale motion the star would be spherically symmetrical in the mean—I am ignoring the possibility of there being a substantial magnetic field—and only the vertical components and of the superadiabatic lapse rate and the convective flux would be nonzero. Moreover, the Reynolds stress would be diagonal; I denote its (3,3) component, the only component that enters hydrostatics when the spherical geometry is ignored, by . In that case the simplest mixing-length estimates for the magnitudes of typical velocity and temperature fluctuations can be obtained by first ignoring the time derivatives and also ignoring the role of the pressure gradient in (7), so that the convective motion is considered to be purely vertical, with velocity . Then one replaces by , where , called the mixing-length, is some characteristic lengthscale of the largest turbulent eddies. For upwardly directed turbulent motion, this leads to where is the thermal diffusivity: . In obtaining these estimates it was assumed that is much less than the scales of variation, , of all the other variables appearing in (13) including . That assumption is the basis of so-called local theories and allows fluctuations at any point to be expressed in terms of the mean state at that point. It is violated in any of the usual mixing-length models of stellar convective envelopes, for which the calibration of , usually to reproduce the solar radius, leads to a value of in excess of unity (e.g., [7, 8]).

To estimate the average one simply multiplies the estimates of and from (13), and introduces a scaling factor to account for the imperfect correlation between velocity and temperature fluctuations, giving where is the product of the Prandtl number and a locally defined Rayleigh number based on the lengthscale . The factor depends on the way in which the mixing-length ideas have been employed to describe the dynamics; in a description of the turbulence more realistic than those assuming that all motion is vertical, it depends also on the supposed geometry of the turbulent eddies. The (3,3) component of the Reynolds stress tensor can be estimated similarly: in which plays a similar role to in (14). The turbulent kinetic-energy flux is usually ignored; it is proportional to , with a constant of proportionality that depends on the (uncertain) spatial correlations of the turbulent velocity field.

Deep in stellar convection zones is very large—in the solar convection zone it exceeds 10^{20}. This implies that the coherent convective motion is essentially adiabatic. Microscopic diffusion of heat is negligible, and
Then
which, of course, are independent of the thermal conductivity .

The formula (18) for has been written as the product of a gradient and a transfer coefficient . The reason for so doing is that it is this form that has been the most common basis for generalization to situations that are no longer spherically symmetrical in the mean. What is sometimes done (e.g., [9]) is to assume that (18) is the 3-component of a vector equation and that the transfer coefficient is a scalar, even though it appears to depend on the vectors and . Thus one obtains with the turbulent heat flux being parallel to the superadiabatic temperature gradient. One might more realistically envisage trying to generalize the formula to something like since generally, as will become apparent later, one would not expect and to be parallel.

Accepting (19) exposes the difficulty of principle in deciding how the coefficient should be computed. Perhaps the simplest hypothesis would be to leave as it stands and consistently replace by the magnitude of in the formula (18) for . But in practice even cruder approximations are often adopted, which at best are based on the value of computed from an initial trial model in which is constrained to be zero, even though the resulting available heat energy advected by the large-scale flow in the final model may turn out to be comparable with .

The treatment of the Reynolds stress in the presence of shear is also motivated by the desire to reduce the formula to a linear Fick law. In most cases the usual scalar viscous law is adopted, with a turbulent shear viscosity derived in some way from the product of and a turbulent-velocity estimate such as that in (17).

It is convenient to consider the Reynolds stress tensor to be divided into two parts: a part that vanishes when and the part that does not. In the absence of a magnetic field or other perturbing agent, is diagonal and acts as a (generally anisotropic) turbulent pressure . It is usually ignored in computations of stellar structure. This may be dangerous, because the presence of any anisotropy in can modify large-scale flow (e.g., [10, 11]).

In modelling the remaining terms in the Reynolds stress, once again an assumption of localness is usually adopted so that it may be presumed that the mean flow may be replaced by a truncated Taylor expansion about any point of interest. Thus one can write There is no term depending on undifferentiated mean velocities because the formula must be invariant under Galilean transformation.

Gross simplifications of the physics have always been made in attempts to evaluate the viscosity tensor. For example, Wasiutyński [12] assumed fluid elements to travel with constant acceleration, relative to an inertial frame, from initial velocities having uncorrelated components. He obtained a formula which in Cartesian coordinates would have been, had spherical geometrical terms been ignored, where is a symmetrical tensor with principal axes in the east-west, north-south, and vertical directions. Such a form, simplified further by assuming axisymmetry about the vertical, was adopted in some of the early discussions of the sun’s equatorial acceleration (e.g., [10, 12, 13]). The assumptions upon which the derivation is based imply, for example, that the convective velocities are not influenced significantly by the solar rotation. Some investigators assume complete isotropy of the turbulent motion, and no resistance to dilatation of the mean flow, which yields where is a scalar. This is the ordinary shear viscosity law, which, unlike the anisotropic formulae, has a (no doubt erroneous) tendency to push the mean flow towards a state of rigid rotation.

#### 4. The Eddy-Ensemble Approach

In order to refine the estimates of the turbulent fluxes one must consider more carefully the dynamics of the turbulent eddies. I adopt here an eddy-ensemble approach developed originally for pulsating stars by Gough [6, 14], who subsequently discussed a generalization to accommodate rotating flows [15] which this paper extends. The model representing the turbulent motion is much more carefully defined than in the previous section, and it is correspondingly more complicated to evaluate. Yet in the simple circumstance of convection in a fluid with no mean background flow it yields essentially the same results. The merit in the approach is that it lays bare the underlying assumptions, thereby making it clear how the theory can be extended to more complicated situations and what additional assumptions, if any, are needed to do so.

It is safer and simpler to compute the fluxes directly from the formulae and rather than to attempt an evaluation of transport coefficients such as and or generalizations of them, as is more common (e.g., [16, 17]). The reason is that (20) and (21) depend on assumptions additional to those of the basic turbulence theory and must therefore predict fluxes that are less soundly based and more difficult to test. Moreover, the tensors and are of higher rank than the fluxes they determine, so one risks exacerbating the complexity of the calculation unnecessarily.

##### 4.1. Axially Symmetric Convection

In order to establish the procedure, I consider first the simple case of a star that is spherically symmetrical. Then the turbulence is (statistically) axially symmetrical about the vertical, and has only two independent components.

The approach I adopt here assumes that an eddy comes into being spontaneously (as a result of the convective instability acting on some random fluctuation in the medium) and subsequently grows in accordance with equations of motion linearized about the mean, background, state, and at all times there is some likelihood that it breaks up due to internal shear instability. The disruption is regarded as being instantaneous and occurs stochastically with probability proportional to the rms rate of strain in the eddy. It is assumed that the smaller-scale motion so generated contributes negligibly to the fluxes. I shall describe here only a local theory, using the Boussinesq approximation with constant within any eddy. The discussion is based on the review by Gough [18], where further elaborations and refinements, such as an approximate account of the smaller-scale motion, can be found.

To render the notation more transparent I rename the coordinates (), use vector notation where convenient, and denote the components of the fluctuating velocity field by (). The convective flow is represented by a superposition of eigenfunctions of the linear problem, each of which is confined between two horizontal planes a distance from its central level . I adopt stress-free isothermal boundary conditions, which yield the simplest formulae. Thus for , a convectively unstable eddy created at time with an initial vertical velocity amplitude may be represented by in which where the planform satisfies and where is the horizontal component, with magnitude , of a wavenumber characterizing the eddy, having magnitude ; is its vertical component. Also , and so forth, is the horizontal Laplacian operator; the overbar denotes horizontal, or ensemble, average. Furthermore, the growth rate of the mode is given by in which and now (e.g., [19]), which differs from the definition (15) by just a constant geometrical factor. Where and (although, of course, there are other eddies there). The formulae can easily be generalized to account for optically thin eddies [6]. Note that formally stable eddies also exist, with , but their direct contribution to the fluxes is not significant; they formally contribute indirectly via a filling factor to be introduced below. The geometrical factor and the planform define the degree of anisotropy of the turbulence.

An eddy is presumed to come into existence as a residue of the turbulent cascade from the breakup of previously existing eddies. The instant at which it is considered to have an existence of its own is not completely specified, although to be unstable it must have a temperature fluctuation that has a (randomly achieved) positive correlation with the vertical velocity. The actual criterion that is adopted to define the initial conditions is absorbed into the definition of the scaling factor appearing in (35) and (37) for the turbulent fluxes, which in practice is subsequently calibrated by comparison with astronomical observation.

The turbulent fluxes can now be computed at height by first constructing the appropriate fluxes due to a single eddy and then averaging over all eddies that intersect the horizontal plane through , weighting the average by the probability that the eddy has not been disrupted. The difficulties in performing this average arise from the computation of that probability, from assigning the distribution of initial amplitudes and the rate at which eddies are created, and for assigning a probability distribution to the eddy geometry characterized by the parameter and, in the case of rotating flow discussed in the next subsection, the function .

Let us consider first the computation of the disruption probability. Eddy disruption is thought to occur predominantly as a result of shear instability within the eddy; its probability of occurring is therefore proportional to the magnitude of that shear. The coefficient of proportionality depends on the geometry of the eddy and is difficult to estimate because shear instability is ill understood. Perhaps the most plausible first guess is to take as the measure of the internal shear the square root of the total squared rate of strain , where and the angular brackets denote average over the volume occupied by the eddy. The ansatz is plausible, for at least it has the correct form in the limits of large and small when the eddy motion approaches rectilinear shear. From the structure (25) of the velocity field, the disruption probability per unit time can then be evaluated to be where is some (constant) parameter to be calibrated later. The probability that an eddy created at time still exists is therefore where . It is an important assumption of the theory that the eddy dynamics is dominated by its growth against the unstably stratified environment, which can be well approximated by linear theory. Then the duration of the process of eddy disruption can be regarded as being short compared with the mean eddy lifetime , which permits disruption to be considered to be instantaneous. This implies that the mean eddy lifetime is given by the final, approximate, expression having been obtained as the leading term in an expansion for small .

We come now to the distribution of and the eddy creation rate. I shall adopt the common practice of assuming that the turbulent convective flow is dominated at any point by eddies with a unique size and shape. Since the flow in a spherically symmetrical (nonrotating) star is statistically axisymmetrical, all orientations of the planform are equally probable. Furthermore, since the theory is local and, for the present, the convection zone is considered to be time independent in the mean, the formula for the creation rate and the initial amplitude depend on neither nor at any given value of .

If is the rate of creation of eddies centred at per unit eddy-volume, then the proportion of the volume occupied by eddies at time is whence ; is essentially the filling factor of eddies with positive linear growth rates.

One is now in a position to evaluate the convective heat flux by integrating the contributions from all eddies. Formally Only eddies that exist at height contribute to the integral over eddy locations , so the outer integral in this simple theory extends only from to . Moreover, consistent with the Boussinesq approximation, when evaluating that integral the spatial variation of can be ignored.

Note that if the horizontal and vertical wavenumbers are of the same order of magnitude, this expression for the probability is similar to that derived from the more literal interpretation of the mixing-length annihilation hypothesis suggested by Spiegel [20] in terms of vertically rising and falling fluid elements, which yields a breakup probability per unit time of , which is proportional to . It is evident from (30) that the two expressions are rather different, especially for the highly elongated eddies having that are sometimes favoured by that literal interpretation [6], particularly if only those eddies with the greatest growth rates are considered, as Spiegel suggested. According to the formula proposed here, such eddies are prone to rapid disintegration by the shear in the highly elongated turbulent flow. Nevertheless, this approach is tantamount to no more than a sophisticated mixing-length formalism. It is evident that the manner in which the mixing-length hypothesis is interpreted affects the predicted anisotropy.

With these assumptions the averages can easily be performed. The analysis is essentially the same as that presented by Gough [6, 18]. The resulting fluxes are which is essentially of the same form as (14). In this formula, which demonstrates how the filling factor , the amplitude at which the eddy is considered to come into existence (which determines the value of ), and the scaling factor defining the eddy destruction rate combine together into a single (calibrateable) parameter. Note that if the assumption that is proportional to a scale height of the background state is made: , the constant of proportionality combines similarly. Equation (35) also demonstrates that the fluxes depend only weakly on , and therefore only weakly on when one considers the eddy to have come into existence. One can evaluate the mean-squared velocities for determining the Reynolds stress similarly: Equations (35)–(37) with (29) are the same as (14)–(16), except for the constants multiplying the functional forms of , , and , and they show the relation between and , and and . Finally one can evaluate the flux of kinetic energy, , likewise: where is a velocity correlation constant. It requires yet another independent assumption to specify . If, for example, it is assumed that and that every horizontal plane is tessellated with Christopherson hexagons (e.g., [21]), which can be represented by six wavenumbers of equal magnitude orientated uniformly in a circle, then (cf. [22, 23]); if instead the plane is filled with closely packed cylinders with stagnant fluid in the interstices, as Böhm-Vitense [2] might prefer, then .

The analysis has not yet provided a method for computing the anisotropy parameter . The simplified formalism presented here is incomplete, and it exposes important issues that must be resolved in order to close what remains a rather naive eddy-ensemble approach. However, before attempting a generalization, an observation of the structure of the formulae (35)–(37) is perhaps not out of place.

Although it has been assumed that turbulent eddies with a given value of dominate the flow, it is appreciated that actually this is merely meant to be representative of a spectrum of eddy shapes. What is meant by an eddy that dominates the flow is one that contributes maximally to a particular flux. Thus is maximized at fixed when , and is maximized at . If these values are used separately for the two fluxes, which is not an implausible procedure if the distribution of eddy shapes is rather flat, one finds This relation is also a property of a single eddy with (which is close to the geometric mean of the two stress-maximizing values, namely ). So perhaps this intermediate value of should be used to describe the representative eddy for computing the stress tensor. It is interesting, though perhaps merely coincidental, that the eddy shape that maximizes at fixed has , which does not differ substantially from 2.

##### 4.2. Convection in a Rotating Fluid

Generalization of the procedure outlined in the previous subsection in the presence of a general background mean flow is not straightforward, because it requires additional hypotheses. Nevertheless, some progress can be made.

The large-scale currents in stars envisaged here include both meridional circulation and zonal flow. I continue to assume that a meaningful separation of the scales of motion is possible, so that a local theory provides a fair first approximation to the turbulent fluxes. As usual only the leading terms in the Taylor expansion of about a point will be considered.

Note first that the mean motion in the vicinity of a point can be decomposed into a translation, a rotation, a pure shearing motion, and an isotropic dilatation. The first does not affect the turbulent motion; the last modifies the turbulent pressure but does not add to it off-diagonal terms; it can be dealt with by generalizing the analysis by Gough [6] for sinusoidally pulsating stars. Rotation modifies the turbulence via the Coriolis forces, and pure shearing motion stretches the eddies; both of these do generate off-diagonal terms. Here I consider the special case where only the local rotation is appreciable; it is the simplest of the steady flows to deal with because there is no shear, and the modifications to the eigenfunctions that the mean flow imposes are relatively straightforward to calculate. Steps towards accounting for rectilinear shearing flow are being taken by Smolec et al. [24, 25].

In order to proceed I shall make two additional simplifying assumptions. Firstly, I assume that surfaces of constant potential temperature (or entropy) are horizontal. This may not be strictly the case, as Durney [9] has pointed out, although the absence of a large pole-equator temperature difference in the sun suggests that at least in that case the assumption may not be a bad first approximation, at least in the upper parts of the convection zone. Secondly, I assume that the vorticity of the mean flow is much less than the eddy growth rate, which (almost) circumvents the issue of how axisymmetry is broken by rotation. This assumption too is easily satisfied in the upper parts of the solar convection zone, but not so well for the slowest eddies that occupy much of the zone beneath. Thus, as with the (sometimes implicit) assumption that the mixing length is much less than the background scale heights, the formulae I derive are not strictly valid, and one should be aware of that when they are used in larger-scale calculations.

As before, the linear eigenfunctions of a plane Boussinesq layer are computed, but now the layer is presumed to rotate with angular velocity . The Cartesian axes are orientated such that lies in the plane. In view of the assumption , there is no question of the rotation stabilizing the convection; the growth rate of the convective mode and its eigenfunction can be expanded in powers of , here up to second order: The leading-order solution is given by which are equivalent to the formula (25) with given by (27)–(29) for . Thermal diffusion is important only in the immediately subphotospheric layers of the star. There the growth rate is high, is very small, and therefore so typically is the influence of the rotation. Deeper down where rotation may matter, the motion is close to being adiabatic. Therefore in calculating the corrections to the structure of the eddies it is expedient to assume adiabaticity, approximating by and omitting some terms in the formulae for the eigenfunctions, in order to maintain a degree of lucidity. It would be quite straightforward to retain the nonadiabatic corrections were the necessity to arise, but I refrain from doing so here because the added complexity of the formulae would render them less easy to appreciate. The rotational modification to the convective mode can accordingly be represented by where is the angle of inclination of to the vertical in the positive direction.

One can now go through the procedure described in Section 4.1 to evaluate the turbulent fluxes. That is quite straightforward except for one issue: the determination of the potential horizontal anisotropy of the turbulent velocity amplitudes resulting from the existence of a preferred direction introduced by the rotation. The assumption that the growth of the eddies satisfies linearized dynamics implies that anisotropy formally arises only from an anisotropy of initial conditions. Since the creation of the unstable convective eddy perturbations is not directly addressed by the formalism, an additional assumption must therefore be made. The assumption that I adopt here is simply that the filling factor is isotropic. When it was introduced in the discussion in Section 4.1 of axisymmetric turbulence I implicitly presumed it to be constant. But now that axisymmetry is lost, that presumption must be relaxed. Here, I adopt the principle of detailed balance, by applying (33) separately in every horizontal direction. That implies that the creation rate is proportional to and therefore shares with the same O() anisotropy. However, it transpires that the anisotropy it imparts to the velocity field is even smaller, reducing the O() terms by an O() factor; it can therefore be neglected. After carrying out the averaging, there results in which and where and are given by (35) and (37). The general formula for is algebraically complicated, but in the case when the planform is composed of randomly orientated harmonic hexagons it simplifies to where is given by (39).

The most obvious influence of the rotation on the convective fluxes is to rotate them principally out of the plane of and , namely, in the direction in the coordinate system adopted here. The angle of rotation is proportional to . There are additional O() components generated in the plane of and , and also the overall magnitudes of the fluxes are reduced somewhat, also by O(). The rotation of the momentum fluxes generates off-diagonal terms in the Reynolds-stress tensor.

#### 5. Discussion

Many assumptions have been necessary to arrive at a tractable procedure. Most notably, it has been assumed here that the turbulent medium can be represented as an ensemble of eddies that not only do not interact directly with each other—they do so implicitly via the background (mean) stratification, however—but that nonlinear advection processes and nonlinearities in the thermodynamics within an eddy (the latter being a consequence of the assumption that the theory is local, which is intimately linked with the Boussinesq assumption) can be ignored throughout almost the entire lifetime of each eddy. The only explicit nonlinearity is the eddy-disruption process, which is assumed to occur on a timescale short compared with the eddy lifetime, which is essential, although not sufficient, for the validity of the linearization assumption: each eddy is presumed to grow via the linear influence of the background state, although there is an implicit nonlinear backreaction on the mean state via the turbulent fluxes. That reaction has been presumed to preserve the alignment of the superadiabatic temperature gradient with gravity, although one can, with additional assumption about the preservation (or otherwise) of eddy shape, take a departure formally into account. Relating that to the fluxes, however, would require a subsequent global calculation.

The formalism is basically a mixing-length approach, which itself is also characterized by linear reasoning, either explicitly when describing the energy exchange with the background state or implicitly via the preservation of eddy size and shape, which is not predicted by the simplest forms of the theory; in any case, it always requires admitting the necessity of an additional assumption about eddy shape. Even though in linear theory the harmonic relation (26) is sufficient for determining the growth rate and the and dependence of the eigenfunctions, more detailed specification is necessary for relating the kinetic-energy flux to the heat flux and the Reynolds stress, in particular for determining the correlation constant . Eddy size is basically the mixing length, which is not an integral part of the theory. Therefore nor is any procedure to define how it changes with changing vorticity (and, in more general circumstances, the changing shear) in the mean flow. Here I have assumed that the “filling factor” of eddies, a fundamental property of the eddy correlation, is unaffected, leaving the way clear to determine the effect on the eddy shape from the linearized dynamics. As always in such a theory, the value of the vertical extent of the eddy is left to other considerations. In particular, no attempt has been made to address how it might change with the properties—in the simple case considered in this paper, the value of the local rotation—of the background flow.

There remains much work to do before an even modestly reliable general theory emerges. Relaxing the restrictions that were imposed in the discussion above not only complicates the analysis but also requires one to face new problems whose resolution may depend on introducing new hypotheses. Accepting that rotation is not small compared with eddy growth rates, for example, raises the issue of the dependence of eddy creation on the orientation of the eddy. Imagining constant entropy surfaces to deviate from the horizontal introduces the possibility of baroclinic instability. How do these compete with the convective modes? It is easy to invent a set of superficially plausible hypotheses for tackling these problems, but some other no less implausible models may yield rather different results. The simple formulae (14)–(16), for example, result from a relatively wide variety of physical pictures of the turbulent motion, but the more subtle results such as (44)–(47) are more sensitive to the details of the assumptions. Can one even hope to develop a theory on the basis of linear eigenmodes? Or do nonlinearities play an essential role during the growth of an eddy? How inaccurate is the process of separating scales? It is certainly quite evident that assuming a unique eddy lengthscale at any location can be a poor recipe, especially near the edges on convection zones where local stability changes sense. One might improve the description of the small-scale turbulence by using a nonlocal theory, but how can one develop a procedure that is both tractable and realistic? Can the mixing-length formalism in any guise be expected to provide an adequate description of the turbulence?

It has been argued that because mixing-length theory is even less reliable at predicting fluxes in complicated situations than it is in the simplest of cases, the complexities of its consequences should be ignored and that it is not worthwhile attempting to generalize the theory. I disagree with this point of view. The mixing-length theory is still the only practical tool available for modelling stellar convection zones for evolutionary calculations—although work such as that by Trampedach and Augustson [26] may change that, but perhaps not before mixing-length theory has been calibrated more precisely [27]—and it should be utilized to fullest advantage until it is superseded. Why prefer a formula such as (23), which is surely incorrect, to one that attempts to include the effects of phenomena that are known to exist? Maybe the results obtained with a more intricate theory will at first be numerically no better at reproducing astronomical observations, but only by studying the consequences of rotation, shear, baroclinicity, and also of magnetic fields on the turbulence can some idea of the importance of these factors be acquired. Ogilvie [28], for example, has made some progress with the latter, using a rather different approach from the one adopted here. The introduction of each new factor tends to require new hypotheses, or to bring into prominence already existing hypotheses upon which the theory did not previously depend in an important way. With these come new parameters which need to be determined. In principle we have at our disposal laboratory experiments and numerical experiments performed with the solar modelling programmes to calibrate the formulae. Both should be utilized—indeed, some have already (e.g., [24, 25, 29]), but only to a limited extent—although eventually care must also be taken to ensure that sufficient checks remain to test the predictive power of the theory.

#### Acknowledgments

The author thanks Pierre Lesaffre for encouraging him to reformulate the theory and extend the expansion to second order in and then to publish it, and he thanks Günter Houdek for useful discussions; also Paula Younger for retyping the original paper. He is grateful to the Leverhulme Trust for an Emeritus Fellowship.

#### References

- E. A. Spiegel and G. Veronis, “On the Boussinesq approximation for a compressible fluid,”
*Astrophysical Journal*, vol. 131, p. 442, 1960,*correction*: vol. 135, p. 665. View at Publisher · View at Google Scholar - E. Böhm-Vitense, “Über die wasserstoffkonvektionszone in sternen verschiedener und effektivtemperaturen und leuchtkräfte,”
*Zeitschrift für Astrophysik*, vol. 46, pp. 108–143, 1958. View at Google Scholar - W. Unno, “Stellar radial pulsation coupled with the convection,”
*Publications of the Astronomical Society of Japan*, vol. 19, p. 140, 1967. View at Google Scholar - J. P. Cox,
*Principles of Stellar Structure*, Gordon & Breach, London, UK, 1967. - R. Kippenhahn and A. Weigert,
*Stellar Structure and Evolution*, Springer, Berlin, Germany, 1990. - D. O. Gough, “Mixing-length theory for pulsating stars,”
*Astrophysical Journal*, vol. 214, no. 1, pp. 196–213, 1977. View at Google Scholar - D. O. Gough and N. O. Weiss, “The calibration of stellar convection theories,”
*Monthly Notices of the Royal Astronomical Society*, vol. 176, pp. 589–607, 1976. View at Google Scholar - J. Christensen-Dalsgaard, W. Däppen, S. V. Ajukov et al., “The current state of solar modeling,”
*Science*, vol. 272, no. 5266, pp. 1286–1292, 1996. View at Google Scholar - B. R. Durney, “On the constancy along cylinders of the angular velocity in the solar convection zone,”
*Astrophysical Journal*, vol. 204, no. 1, pp. 589–596, 1976. View at Publisher · View at Google Scholar - R. Kippenhahn, “Differential rotation in stars with convective envelopes,”
*Astrophysical Journal*, vol. 137, p. 664, 1963. View at Publisher · View at Google Scholar - E. A. Spiegel and J.-P. Zahn, “The solar tachocline,”
*Astronomy and Astrophysics*, vol. 265, no. 1, pp. 106–114, 1992. View at Google Scholar - J. Wasiutyński, “Studies in hydrodynamics and structure of stars and planets,”
*Astrophisica Norvegica*, vol. 4, 1946. View at Google Scholar - W. J. Cocke, “On the solar differential rotation and meridional currents,”
*Astrophysical Journal*, vol. 150, p. 1041, 1967. View at Publisher · View at Google Scholar - D. O. Gough, “The pulsational stability of a convective atmosphere,” in
*Geophysical Fluid Dynamics II*, pp. 49–84, Woods Hole Oceanographic Institution, 1965. View at Google Scholar - D. O. Gough, “The estimation of fluxes due to small-scale turbulence,” in
*Proceedings of the Workshop on Solar Rotation*, G. Belvedere and L. Paternò, Eds., pp. 337–350, Università di Catania, Catania, Italy, 1978. - B. R. Durney and H. C. Spruit, “On the dynamics of stellar convection zones —the effect of rotation on the turbulent viscosity and conductivity,”
*Astrophysical Journal*, vol. 234, pp. 1067–1078, 1979. View at Publisher · View at Google Scholar - G. Rüdiger, P. Egorov, L. L. Kitchatinov, and M. Küker, “The eddy heat-flux in rotating turbulent convection,”
*Astronomy and Astrophysics*, vol. 431, no. 1, pp. 345–352, 2005. View at Publisher · View at Google Scholar - D. O. Gough, “The current state of stellar mixing-length theory,” in
*Proceedings of the International Astronomical Union Colloquium*, E. A. Spiegel and J.-P. Zahn, Eds., pp. 15–56, Springer, Heidelberg, Germany, 1977. - P. Ledoux, M. Schwarzschild, and E. A. Spiegel, “On the spectrum of trurbulent convection,”
*Astrophysical Journal*, vol. 133, p. 184, 1961. View at Publisher · View at Google Scholar - E. A. Spiegel, “A generalization of the mixing-length theory of turbulent convection,”
*Astrophysical Journal*, vol. 138, p. 216, 1963. View at Publisher · View at Google Scholar - S. Chandrasekhar,
*Hydrodynamic and Hydromagnetic Stability*, Oxford University Press, Oxford, UK, 1961. - P. H. Roberts, “On non-linear Bénard convection (with an Appendix by K. Stewartson),” in
*Non-Equilibrium Thermodynamics, Variational Techniques, and Stability*, R. Donnelly, R. Herman, and I. Prigogine, Eds., p. 125, University of Chicago Press, 1966. - D. O. Gough, E. A. Spiegel, and J. Toomre, “Modal equations for cellular convection,”
*Journal of Fluid Mechanics*, vol. 68, no. 4, pp. 695–719, 1975. View at Google Scholar - R. Smolec, G. Houdek, and D. O. Gough, “Modelling turbulent fluxes due to thermal convection in rectilinear shearing flow,” in
*Proceedings of the International Astronomical Union Symposium*, N. H. Brummell, A. S. Brun, M. S. Miesch, and Y. Ponty, Eds., vol. 271, pp. 397–398, Cambridge University Press, 2011. View at Publisher · View at Google Scholar - R. Smolec, G. Houdek, and D. O. Gough, “Nonlocal model for the turbulent fluxes due to thermal convection in rectilinear shearing flow,” in
*Proceedings of the 20th Stellar Pulsation Conference Series: Impact of New Instrumentation & New Insights in Stellar Pulsations*, Granada, Spain, 2011. - R. Trampedach and K. Augustson, “Using simulations of solar surface convection as boundary conditions on global simulations,” in
*Proceedings of the International Astronomical Union Symposium*, N. H. Brummell, A. S. Brun, M. S. Miesch, and Y. Ponty, Eds., vol. 271, pp. 403–404, Cambridge University Press, 2011. View at Publisher · View at Google Scholar - R. Trampedach and R. F. Stein, “The mass mixing length in convective stellar envelopes,”
*Astrophysical Journal*, vol. 731, no. 2, article 78, 2011. View at Publisher · View at Google Scholar - G. I. Ogilvie, “On the dynamics of magnetorotational turbulent stresses,”
*Monthly Notices of the Royal Astronomical Society*, vol. 340, no. 3, pp. 969–982, 2003. View at Publisher · View at Google Scholar - P. Tooth and D. O. Gough, “Towards an independent calibration of the mixing-length theory,” in
*Seismology of the Sun and the Sun-Like Stars*, E. J. Rolfe, Ed., ESA, Noordwijk, The Netherlands, 1988. View at Google Scholar