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

Authors: D. O. Gough
Journal: ISRN Astronomy and Astrophysics
 (2012)

---

## 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.


---

## Body

## 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 fieldsv and T via the application of an averaging procedure, which I denote by an overbar. Thus(1)v=v¯+u,T=T¯+ϑ.
I shall refer to fields such as v¯ and T¯ as mean quantities; the fields u 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 coordinatesxi, with x3 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,(2)∂∂t(ρvi)+∂∂xj(ρvivj)=-∂p∂xi-gρδ3i,(3)∂∂t(ρϵ+12ρvivi)+∂∂xi(Fi+ρhvi+12ρvjvjvi)=-gρv3,(4)∂ρ∂t+∂∂xi(ρvi)=0,
where ρ and p are density and pressure, ϵ and h are specific internal energy and enthalpy (viz., internal energy and enthalpy per unit mass), g is the magnitude of the gravitational acceleration gi, Fi is the radiative heat flux, t is time, and δij 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 vi to relate the rate of change of the kinetic energy to the rate of working, with the thermal-energy equation(5)ρcv(∂T∂t+vi∂T∂xi)+p∂vi∂xi+∂Fi∂xi=0,
which is a statement of the first law of thermodynamics. In (5), cv is the specific heat at constant volume. These equations must, of course, be supplemented by an equation of state and an equation determining Fi. Viscous terms have been ignored; so too have perturbations to g.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 v̅i=0, 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(6)∂∂t(ρ̅v̅i)+∂∂xj(ρ̅v̅iv̅j)+∂Rij∂xj=-∂p̅∂xi-gρ̅δ3i,∂∂t(ρ̅ϵ̅+12ρ̅v̅jv̅j)+∂∂xi(F¯i+ρ̅h̅v̅i+Fci+12ρ̅(v̅jv̅j+ujuj¯)v̅i+Kti)=-gρ̅v̅3,∂ρ̅∂t+∂∂xi(ρ̅v̅i)=0.
The only fluxes arising from the small-scale turbulent motion that need to be calculated for closing these equations are evidently the Reynolds stress Rij=ρ̅uiuj¯, the turbulent heat flux Fci=ρ̅h′ui¯, where the prime denotes a turbulent fluctuation, and the turbulent kinetic-energy flux Kti=(1/2)ρ̅ujujui¯.In the absence of a large-scale mean flow the equations for the fluctuations are given in Boussinesq approximation by [1](7)ρ̅(∂ui∂t+uj∂ui∂xj)=-∂p′∂xi+gρ̅δ̅T̅ϑδi3,(8)∂ϑ∂t+ui∂ϑ∂xi-βiui=-1ρ̅c̅p∂Fi′∂xi,(9)∂ui∂xi=0,
where cp is the specific heat at constant pressure and(10)βi≡δ̅ρ̅c̅p∂p̅∂xi-∂T̅∂xi
is the superadiabatic lapse rate; δ=-(∂lnp/∂lnρ)T 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,(11)Fi′=-K∂ϑ∂xi,
where K=4acT̅3/3χ̅ρ̅, in which a is the first radiation constant, c 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(12)Fci=ρcpϑui¯.
(It is perhaps worth pointing out that for a perfect unionizing gas the pressure fluctuation makes precisely no contribution to h′ 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 Fc 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 pt. 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 w. Then one replaces ∂/∂xi by ℓ-1, where ℓ, called the mixing-length, is some characteristic lengthscale of the largest turbulent eddies. For upwardly directed turbulent motion, this leads to(13)w2≃gδlTϑ,(l-1ϑ-β)w≃-κl-2ϑ,
where κ is the thermal diffusivity: κ=K/ρcp. In obtaining these estimates it was assumed that ℓ is much less than the scales of variation, H, 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 ℓ/H in excess of unity (e.g., [7, 8]).To estimate the averagewϑ¯ one simply multiplies the estimates of w and ϑ from (13), and introduces a scaling factor αH to account for the imperfect correlation between velocity and temperature fluctuations, giving(14)Fc≃18αHS-1(4S+1-1)3Kβ,
where(15)S=gδT-1βl4κ2
is the product of the Prandtl number and a locally defined Rayleigh number based on the lengthscale ℓ. The factor αH 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 pt of the Reynolds stress tensor Rij can be estimated similarly:(16)pt≃αRρw2≃14αR(4S+1-1)2ρ(κl)2=αR(1+14S-1-12S-1/2)2ρgδβTl2,
in which αR plays a similar role to αH in (14). The turbulent kinetic-energy flux is usually ignored; it is proportional to ρw3, with a constant of proportionality that depends on the (uncertain) spatial correlations of the turbulent velocity field.Deep in stellar convection zonesS is very large—in the solar convection zone it exceeds 1020. This implies that the coherent convective motion is essentially adiabatic. Microscopic diffusion of heat is negligible, and(17)w≃gδTlβ1/2,ϑ≃lβ.
Then(18)Fc≃αHρcpgδTl2β3/2=αHS1/2Kβ=:Kcβ,pt≃αRgδρTl2β,
which, of course, are independent of the thermal conductivity K.The formula (18) for Fc has been written as the product of a gradient β and a transfer coefficient Kc. 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 Kc is a scalar, even though it appears to depend on the vectors βi and gi. Thus one obtains(19)Fci=Kcβi,
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(20)Fci=Kcijβj,
since generally, as will become apparent later, one would not expect Fci and βi to be parallel.Accepting (19) exposes the difficulty of principle in deciding how the coefficient Kc should be computed. Perhaps the simplest hypothesis would be to leave g as it stands and consistently replace β by the magnitude (βiβi)1/2 of βi in the formula (18) for Kc. But in practice even cruder approximations are often adopted, which at best are based on the value of Kc computed from an initial trial model in which v̅i 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 Fci.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 tensorRij to be divided into two parts: a part R1ij that vanishes when v̅i=0 and the part R0ij that does not. In the absence of a magnetic field or other perturbing agent, R0ij is diagonal and acts as a (generally anisotropic) turbulent pressure pti. It is usually ignored in computations of stellar structure. This may be dangerous, because the presence of any anisotropy in pti 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(21)ρuiuj¯=pt(i)δij+Mijkm∂v̅m∂xk+higher-order terms.
There is no term depending on undifferentiated mean velocities v̅i 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,(22)Mijkm=(δinδjm+δimδjn)Mnk,
where Mnk 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(23)Mijkm=(δikδjm+δimδjk+δijδkm)η,
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ρuiuj¯ and ρcpϑui¯ rather than to attempt an evaluation of transport coefficients such as Kcij and Mijkm 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 Kcij and Mijkm 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, andRij=pt(i)δij 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βi 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 thexi coordinates (x,y,z), use vector notation where convenient, and denote the components ui of the fluctuating velocity fieldu by (u,v,w). 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 ±(1/2)ℓ from its central level z0. I adopt stress-free isothermal boundary conditions, which yield the simplest formulae. Thus for z0-(1/2)ℓ≤z≤z0+(1/2)ℓ, a convectively unstable eddy created at time t0 with an initial vertical velocity amplitude W0 may be represented by(24)(u,v,w,ϑ)=(û,v̂,ŵ,ϑ̂)W0eq(t-t0),
in which(25)û(x,y,z-z0)=-kh-2kvsin⁡[kv(z-z0)]fx,v̂(x,y,z-z0)=-kh-2kvsin⁡[kv(z-z0)]fy,ŵ(x,y,z-z0)=cos⁡[kv(z-z0)]f,ϑ̂(x,y,z-z0)=βq-1ŵ(x,y,z-z0),
where the planform f(x,y) satisfies(26)∇h2f+kh2f=0,f2¯=1,
and where kh=(kx,ky,0) is the horizontal component, with magnitude kh, of a wavenumber k characterizing the eddy, having magnitude k; kv=π/ℓ is its vertical component. Also fx:=∂f/∂x, and so forth, ∇h2 is the horizontal Laplacian operator; the overbar denotes horizontal, or ensemble, average. Furthermore, the growth rate q of the mode is given by(27)q=-12k2κ+μ2+14k4κ2=12(4S+1-1)k2κ,
in which(28)μ2=gδβΦT,Φ=1+kv2kh2,
and now(29)S=μ2κ2k4
(e.g., [19]), which differs from the definition (15) by just a constant geometrical factor. Where |z-z0|>(1/2)ℓ,u=0and ϑ=0 (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 ϑ̂=-βq-1ŵ, but their direct contribution to the fluxes is not significant; they formally contribute indirectly via a filling factor f̂ to be introduced below. The geometrical factor Φ and the planform f 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 heightz by first constructing the appropriate fluxes due to a single eddy and then averaging over all eddies that intersect the horizontal plane through z, 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 f.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 roote=e0exp[q(t-t0)] of the total squared rate of strain e2=〈eijeij〉, where eij=(1/2)(∂ui/∂xj+∂uj/∂xi) 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 kv/kh 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(30)γe=12γ(1-Φ-1)kvW0exp⁡[q(t-t0)],
where γ is some (constant) parameter to be calibrated later. The probability that an eddy created at time t0 still exists is therefore(31)P(t;t0)=exp⁡[-λ(eq(t-t0)-1)],
where λ=γq-1e0. 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(32)τ=γe0∫t0∞(t-t0)Pdt≃ln⁡(1λ)q-1,
the final, approximate, expression having been obtained as the leading term in an expansion for small λ.We come now to the distribution ofW0 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 f 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 W0 depend on neither z0 nor t0 at any given value of z.Ifn(z0) is the rate of creation of eddies centred at z=z0 per unit eddy-volume, then the proportion of the volume occupied by eddies at time t is(33)f̂=n∫-∞tPdt0=nτ,
whence n=f̂τ-1; f̂ 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(34)Fc=ρcp∫n(∫-∞twϑPdt0)dz0.
Only eddies that exist at height z contribute to the integral over eddy locations z0, so the outer integral in this simple theory extends only from z-ℓ/2 to z+ℓ/2. Moreover, consistent with the Boussinesq approximation, when evaluating that integral the spatial variation of n 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 w/ℓ, which is proportional to kvW0. It is evident from (30) that the two expressions are rather different, especially for the highly elongated eddies having kv≪kh 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(35)Fc=ΛΦ-5/2(Φ-1)ρcpTgδq3l2=ΛΦ-5/2(Φ-1)(1+4S-1-12S-1)3ρcpgδβTl2β,
which is essentially of the same form as (14). In this formula,(36)Λ=4f̂π2γ2ln⁡1/λ,
which demonstrates how the filling factor f̂, the amplitude at which the eddy is considered to come into existence (which determines the value of e0), 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 H of the background state is made: ℓ=αH, 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:
(37)pt=ρw2¯=ΛΦ-3(Φ-1)ρl2q2=ΛΦ-3(Φ-1)(1+4S-1-12S-1)2ρgδβTl2,(38)u2¯=v2¯=12(Φ-1)w2¯.
Equations (35)–(37) with (29) are the same as (14)–(16), except for the constants multiplying the functional forms of Fc, pt, and S, and they show the relation between αH and αR, and Λ and Φ. Finally one can evaluate the flux of kinetic energy, Kt, likewise:(39)Kt=ΛCΦ-9/2(Φ-1)3/2(Φ+3)(1+4S-1-12S-1)3×ρ(gδβT)3/2l3,
where C=82f3¯/3π2γ is a velocity correlation constant. It requires yet another independent assumption to specify C. If, for example, it is assumed that γ=1 and that every horizontal plane is tessellated with Christopherson hexagons (e.g., [21]), which can be represented by six wavenumbers kh of equal magnitude k orientated uniformly in a circle, then C=4/33π2≃0.078 (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 C≃0.034.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 w2¯ is maximized at fixed ℓ when Φ=1.5, and u2¯ is maximized at Φ=3. 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(40)u2¯=12w2¯.
This relation is also a property of a single eddy with Φ=2 (which is close to the geometric mean of the two stress-maximizing values, namely 3/2=2.1). 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 Fc at fixed ℓ has Φ=5/3, 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 ofv¯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Ω=(1/2)curlv¯. The Cartesian axes are orientated such that Ωlies in the x-z plane. In view of the assumption Ω:=|Ω|≪q≃μ, 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:(41)q=q0(1+εq̃1+ε2q̃2+⋯),û=û0+εû1+ε2û2+⋯,etc.
The leading-order solution is given by(42)ŵ0=cos⁡[kv(z-z0)]f=q0β-1ϑ̂0,û0=-kh-2kvsin⁡[kv(z-z0)]fx=fxfy-1v̂0,
which are equivalent to the formula (25) with q0 given by (27)–(29) for q. 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 q0 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(43)q̃1=0,ŵ1=0,ϑ̂1=0,û1=-2kh-2[cos⁡θsin⁡[kv(z-z0)]fy-sin⁡θcos⁡[kv(z-z0)]fxy],v̂1=2kh-2[cos⁡θsin⁡[kv(z-z0)]fx-sin⁡θcos⁡[kv(z-z0)]fxx],q̃2=-2k-2(fx2¯sin⁡θ+kv2co⁡s2θ)ŵ2=-2sin⁡2θ⋅(z-z0)cos⁡[kv(z-z0)]fx,ϑ̂2=μ-1β[k-2(-sin⁡2θfxx+kv2cos⁡2θf)-ŵ2],û2=-2kh-2sin⁡2θ{cos⁡[kv(z-z0)]-kv(z-z0)sin⁡[kv(z-z0)]}fxx=fxxfxy-1v̂2,
where θ is the angle of inclination of Ωto the vertical in the positive x direction.One can now go through the procedure described in Section4.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 f̂ 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 n is proportional to q and therefore shares with q the same O(ε2) anisotropy. However, it transpires that the anisotropy it imparts to the velocity field is even smaller, reducing the O(ε2) terms by an O(ln1/ε) factor; it can therefore be neglected. After carrying out the averaging, there results(44)Fc=(ε2sin⁡2θ,εsin⁡θ,1)ψFc,(45)Rij=(12(Φ-1)+12ε2Rx02ε2sin⁡2θ012(Φ-1)+12ε2Ryεsin⁡θ2ε2sin⁡2θεsin⁡θ1)ψpt,
in which(46)Rx=Φ-1,Ry=(Φ-1)cos⁡2θ+34sin⁡2θ,ψ=1-2εΦ-1sin⁡θ-ε2Φ-2[6(Φ-1)cos⁡2θ-sin⁡2θ],
and where Fc and pt are given by (35) and (37). The general formula for Kt is algebraically complicated, but in the case when the planform is composed of randomly orientated harmonic hexagons it simplifies to(47)Kt=136π(112ε2(Φ-1)sin⁡2θ,4(Φ-1)εsin⁡θ,1)ψKt,
where Kt is given by (39).The most obvious influence of the rotation on the convective fluxes is to rotate them principally out of the plane ofg and Ω, namely, in the y direction in the coordinate system adopted here. The angle of rotation is proportional to ε=Ω/μ. There are additional O(ε2) components generated in the plane of g and Ω, and also the overall magnitudes of the fluxes are reduced somewhat, also by O(ε2). The rotation of the momentum fluxes generates off-diagonal terms in the Reynolds-stress tensor.


## 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, andRij=pt(i)δij 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βi 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 thexi coordinates (x,y,z), use vector notation where convenient, and denote the components ui of the fluctuating velocity fieldu by (u,v,w). 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 ±(1/2)ℓ from its central level z0. I adopt stress-free isothermal boundary conditions, which yield the simplest formulae. Thus for z0-(1/2)ℓ≤z≤z0+(1/2)ℓ, a convectively unstable eddy created at time t0 with an initial vertical velocity amplitude W0 may be represented by(24)(u,v,w,ϑ)=(û,v̂,ŵ,ϑ̂)W0eq(t-t0),
in which(25)û(x,y,z-z0)=-kh-2kvsin⁡[kv(z-z0)]fx,v̂(x,y,z-z0)=-kh-2kvsin⁡[kv(z-z0)]fy,ŵ(x,y,z-z0)=cos⁡[kv(z-z0)]f,ϑ̂(x,y,z-z0)=βq-1ŵ(x,y,z-z0),
where the planform f(x,y) satisfies(26)∇h2f+kh2f=0,f2¯=1,
and where kh=(kx,ky,0) is the horizontal component, with magnitude kh, of a wavenumber k characterizing the eddy, having magnitude k; kv=π/ℓ is its vertical component. Also fx:=∂f/∂x, and so forth, ∇h2 is the horizontal Laplacian operator; the overbar denotes horizontal, or ensemble, average. Furthermore, the growth rate q of the mode is given by(27)q=-12k2κ+μ2+14k4κ2=12(4S+1-1)k2κ,
in which(28)μ2=gδβΦT,Φ=1+kv2kh2,
and now(29)S=μ2κ2k4
(e.g., [19]), which differs from the definition (15) by just a constant geometrical factor. Where |z-z0|>(1/2)ℓ,u=0and ϑ=0 (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 ϑ̂=-βq-1ŵ, but their direct contribution to the fluxes is not significant; they formally contribute indirectly via a filling factor f̂ to be introduced below. The geometrical factor Φ and the planform f 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 heightz by first constructing the appropriate fluxes due to a single eddy and then averaging over all eddies that intersect the horizontal plane through z, 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 f.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 roote=e0exp[q(t-t0)] of the total squared rate of strain e2=〈eijeij〉, where eij=(1/2)(∂ui/∂xj+∂uj/∂xi) 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 kv/kh 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(30)γe=12γ(1-Φ-1)kvW0exp⁡[q(t-t0)],
where γ is some (constant) parameter to be calibrated later. The probability that an eddy created at time t0 still exists is therefore(31)P(t;t0)=exp⁡[-λ(eq(t-t0)-1)],
where λ=γq-1e0. 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(32)τ=γe0∫t0∞(t-t0)Pdt≃ln⁡(1λ)q-1,
the final, approximate, expression having been obtained as the leading term in an expansion for small λ.We come now to the distribution ofW0 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 f 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 W0 depend on neither z0 nor t0 at any given value of z.Ifn(z0) is the rate of creation of eddies centred at z=z0 per unit eddy-volume, then the proportion of the volume occupied by eddies at time t is(33)f̂=n∫-∞tPdt0=nτ,
whence n=f̂τ-1; f̂ 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(34)Fc=ρcp∫n(∫-∞twϑPdt0)dz0.
Only eddies that exist at height z contribute to the integral over eddy locations z0, so the outer integral in this simple theory extends only from z-ℓ/2 to z+ℓ/2. Moreover, consistent with the Boussinesq approximation, when evaluating that integral the spatial variation of n 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 w/ℓ, which is proportional to kvW0. It is evident from (30) that the two expressions are rather different, especially for the highly elongated eddies having kv≪kh 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(35)Fc=ΛΦ-5/2(Φ-1)ρcpTgδq3l2=ΛΦ-5/2(Φ-1)(1+4S-1-12S-1)3ρcpgδβTl2β,
which is essentially of the same form as (14). In this formula,(36)Λ=4f̂π2γ2ln⁡1/λ,
which demonstrates how the filling factor f̂, the amplitude at which the eddy is considered to come into existence (which determines the value of e0), 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 H of the background state is made: ℓ=αH, 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:
(37)pt=ρw2¯=ΛΦ-3(Φ-1)ρl2q2=ΛΦ-3(Φ-1)(1+4S-1-12S-1)2ρgδβTl2,(38)u2¯=v2¯=12(Φ-1)w2¯.
Equations (35)–(37) with (29) are the same as (14)–(16), except for the constants multiplying the functional forms of Fc, pt, and S, and they show the relation between αH and αR, and Λ and Φ. Finally one can evaluate the flux of kinetic energy, Kt, likewise:(39)Kt=ΛCΦ-9/2(Φ-1)3/2(Φ+3)(1+4S-1-12S-1)3×ρ(gδβT)3/2l3,
where C=82f3¯/3π2γ is a velocity correlation constant. It requires yet another independent assumption to specify C. If, for example, it is assumed that γ=1 and that every horizontal plane is tessellated with Christopherson hexagons (e.g., [21]), which can be represented by six wavenumbers kh of equal magnitude k orientated uniformly in a circle, then C=4/33π2≃0.078 (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 C≃0.034.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 w2¯ is maximized at fixed ℓ when Φ=1.5, and u2¯ is maximized at Φ=3. 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(40)u2¯=12w2¯.
This relation is also a property of a single eddy with Φ=2 (which is close to the geometric mean of the two stress-maximizing values, namely 3/2=2.1). 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 Fc at fixed ℓ has Φ=5/3, 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 ofv¯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Ω=(1/2)curlv¯. The Cartesian axes are orientated such that Ωlies in the x-z plane. In view of the assumption Ω:=|Ω|≪q≃μ, 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:(41)q=q0(1+εq̃1+ε2q̃2+⋯),û=û0+εû1+ε2û2+⋯,etc.
The leading-order solution is given by(42)ŵ0=cos⁡[kv(z-z0)]f=q0β-1ϑ̂0,û0=-kh-2kvsin⁡[kv(z-z0)]fx=fxfy-1v̂0,
which are equivalent to the formula (25) with q0 given by (27)–(29) for q. 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 q0 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(43)q̃1=0,ŵ1=0,ϑ̂1=0,û1=-2kh-2[cos⁡θsin⁡[kv(z-z0)]fy-sin⁡θcos⁡[kv(z-z0)]fxy],v̂1=2kh-2[cos⁡θsin⁡[kv(z-z0)]fx-sin⁡θcos⁡[kv(z-z0)]fxx],q̃2=-2k-2(fx2¯sin⁡θ+kv2co⁡s2θ)ŵ2=-2sin⁡2θ⋅(z-z0)cos⁡[kv(z-z0)]fx,ϑ̂2=μ-1β[k-2(-sin⁡2θfxx+kv2cos⁡2θf)-ŵ2],û2=-2kh-2sin⁡2θ{cos⁡[kv(z-z0)]-kv(z-z0)sin⁡[kv(z-z0)]}fxx=fxxfxy-1v̂2,
where θ is the angle of inclination of Ωto the vertical in the positive x direction.One can now go through the procedure described in Section4.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 f̂ 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 n is proportional to q and therefore shares with q the same O(ε2) anisotropy. However, it transpires that the anisotropy it imparts to the velocity field is even smaller, reducing the O(ε2) terms by an O(ln1/ε) factor; it can therefore be neglected. After carrying out the averaging, there results(44)Fc=(ε2sin⁡2θ,εsin⁡θ,1)ψFc,(45)Rij=(12(Φ-1)+12ε2Rx02ε2sin⁡2θ012(Φ-1)+12ε2Ryεsin⁡θ2ε2sin⁡2θεsin⁡θ1)ψpt,
in which(46)Rx=Φ-1,Ry=(Φ-1)cos⁡2θ+34sin⁡2θ,ψ=1-2εΦ-1sin⁡θ-ε2Φ-2[6(Φ-1)cos⁡2θ-sin⁡2θ],
and where Fc and pt are given by (35) and (37). The general formula for Kt is algebraically complicated, but in the case when the planform is composed of randomly orientated harmonic hexagons it simplifies to(47)Kt=136π(112ε2(Φ-1)sin⁡2θ,4(Φ-1)εsin⁡θ,1)ψKt,
where Kt is given by (39).The most obvious influence of the rotation on the convective fluxes is to rotate them principally out of the plane ofg and Ω, namely, in the y direction in the coordinate system adopted here. The angle of rotation is proportional to ε=Ω/μ. There are additional O(ε2) components generated in the plane of g and Ω, and also the overall magnitudes of the fluxes are reduced somewhat, also by O(ε2). 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 q and the z and t 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 C. 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.



---

*Source: 987275-2012-02-29.xml*