text:
# Closed Trapped Surfaces in Cosmology

## 1 Introduction

Since Roger Penrose’ pioneering paper of 1965 , the existence of closed trapped surfaces (‘CTSs’) has been understood as a geometrical condition that, jointly with suitable energy conditions, in many circumstances leads to the existence of space-time singularities in the context of both gravitational collapse and cosmology. This understanding has been codified in the series of singularity theorems proved by Penrose and Stephen Hawking . The various theorems involve different combinations of geometric requirements and energy conditions. Under the assumptions of standard hot big bang theory, these conditions will indeed be met in the cosmological context, because the existence of the black-body cosmic background radiation implies the existence of CTSs in the era between decoupling and the present day  and so leads to prediction of a (classical) singularity at the start of the universe. However it is now known that scalar fields can violate some of the energy conditions, thereby providing the foundation of the inflationary universe paradigm . The possibility then arises of singularity avoidance in realistic early universe models because of these energy condition violations, and indeed even avoidance of a quantum gravity regime is possible , despite the existence of CTSs. The purpose of this paper is to revisit the relation between CTSs and spacetime singularities in this cosmological context, characterising cases where existence of CTSs do not imply the existence of singularities. In particular we examine the case of the de Sitter Universe, showing that CTSs exist in these spacetimes even though (in maximally extended form) they are both geodesically complete and stable to perturbations.

### 1.1  Closed Trapped Surfaces

Consider a spacelike 2-surface S with spherical topology. It is a  past closed trapped surface if both families of past null geodesics orthogonal to the 2-surface S are converging, i.e. if their divergences (evaluated at the surface) are negative. Similarly future closed trapped surfaces occur if we replace ”past” by ”future” in the above. Both past and future closed trapped surfaces will be referred to as  closed trapped surfaces (‘CTSs’). We will in this paper be concerned with 2-surfaces that are 2-spheres with a group of isometries transitive on them, so they are homogeneous 2-dimensional subspaces of spacetime. Then the value of the divergence of each family of normals is constant over the 2-surface, and can be characterised by a single number on each 2-sphere. Thus we will in fact be considering existence of homogeneous closed trapped 2-spheres.

Marginally closed trapped surfaces exist if the divergences are non-positive for these families of null geodesics; that is, if the divergences are either zero or negative rather than strictly negative.

### 1.2 Energy Conditions

Energy conditions conditions generically lead to convergences of irrotational familes of non-spacelike and null geodesics respectively.

#### 1.2.1 Non-spacelike convergence condition

This is the condition

$$R_{ab}K^aK^b0\mathrm{for}\mathrm{all}\mathrm{non}\mathrm{spacelike}\mathrm{vectors}K^a.$$

For perfect fluids, this translates into $`\mu +p0,`$ $`\mu +3p0,`$ which will be true for all ordinary matter. For scalar fields, it becomes $`\frac{1}{2}\dot{\varphi }^20,\dot{\varphi }^2V(\varphi ),`$ hence is violated when the slow rolling condition $`\dot{\varphi }^2V(\varphi )`$ is satisfied. A cosmological constant is the case $`\dot{\varphi }^2=0,V>0`$ and hence violates this condition.

#### 1.2.2 Null convergence condition

This is the condition

$$R_{ab}K^aK^b0\mathrm{for}\mathrm{all}\mathrm{null}\mathrm{vectors}K^a$$
(1)
which is implied by the previous:

$$\{\mathrm{Non}\mathrm{spacelike}\mathrm{convergence}\mathrm{condition}\}\{\mathrm{Null}\mathrm{convergence}\mathrm{condition}\}.$$

For perfect fluids, this translates into $`\mu +p0,`$ which will be true for all ordinary matter. For scalar fields, it is $`\frac{1}{2}\dot{\varphi }^20,`$ and so is true for all ordinary scalar fields (we discount the possibility of ‘phantom matter’ that violates this condition, see Gibbons  for a discussion).

### 1.3 Focussing

The equation determining the evolution of the convergence $`\theta =K_{;a}^a`$ of hypersurface-orthogonal null geodesics is

$$\frac{d\theta }{dv}+\frac{1}{2}\theta ^2=R_{ab}K^aK^b2\sigma ^2,$$
(2)
while the shear propagation equation is

$$\frac{d}{dv}\sigma _{mn}=\theta \sigma _{mn}C_{manb}K^aK^b.$$

This shows that the shear can only remain zero either for very special spacetimes (e.g. Robertson-Walker spacetimes where $`C_{manb}=0)`$, or for very special null rays in a more generic spacetime, so that $`C_{manb}K^aK^b=0`$ at every point on the null geodesics because the geodesic tangent vector is in a special relation to the Weyl tensor (it is a principal null direction). Once the shear is non-zero, it acts as a source term in the null Raychaudhuri equation (2).

When the null convergence condition (1) is true, there is an exceptional case and a generic case for familes of hypersurface orthogonal null geodesics. The exceptional case occurs if $`\theta =0,`$i.e. no focussing occurs:

$$\theta =0R_{ab}K^aK^b=0,\sigma ^2=0C_{acbd}K^aK^b=0.$$

This cannot be true in a generic cosmological context, for example a perturbed Robertson-Walker universe, when both the Ricci and Weyl tensor conditions will be violated along a generic null ray. The generic case is when either $`\theta _0<0`$ or $`\theta _0>0`$. Both imply $`\theta \mathrm{}`$ within a finite affine distance, either to the past or the future. Then the null rays intersect at a caustic, so the surface generated by the null geodesics experiences self-intersections either before or at those events.

Hence when we consider CTSs in realistic cosmologies, which will always satisfy the null convergence condition, both families of converging orthogonal null geodesics, which generate the boundary of their past for the case of past CTSs, will self-intersect. Then by well-known causal theorems, these null rays will lie from then on inside the pasts of the CTSs; the boundary of the past is therefore compact. This is what underlies the singularity theorems.

### 1.4 The Major Singularity Theorems

The major singularity theorems referring to CTSs are given in Hawking and Ellis  (‘HE’). Each case assumes a CTS or roughly equivalent condition, plus the following:

HE Theorem 1  \- A non-compact Cauchy surface and the null convergence condition,

HE Theorem 2  \- The non-spacelike convergence condition and a causality condition,

HE Theorem 3  \- The non-spacelike convergence condition and a causality condition.

Roughly speaking: in each case the boundary of the past comes to an end because of existence of self-intersections points in its generating geodesics; but the past is contained within the boundary, hence a singularity must occur.

We now consider how these theorems apply to Friedmann-Lemaitre (‘FL’) model universes, interpreted here as universes with a Robertson-Walker geometry and matter content of ordinary matter and/or a scalar field.

## 2 Friedmann-Lemâitre Models

Friedmann-Lemaitre universe have a Robertson-Walker (‘RW’) metric which can be represented in the form

$$ds^2=dt^2+S^2(t)\left(dr^2+f^2(r)(d\theta ^2+\mathrm{sin}^2\theta d\varphi ^2)\right)$$
(3)
The scale factor is $`S(t)`$, the matter 4-velocity is: $`u^a=\delta _0^au_a=g_{ab}u^b=\delta _a^0,`$ and $`f(r)=(\mathrm{sin}r,r,\mathrm{sinh}r)`$ if ($`k=1,0,1)`$respectively.

The metric determinant $`g`$ is

$$g=S^6(t)f^4(r)\mathrm{sin}^2\theta \sqrt{g}=S^3(t)f^2(r)\mathrm{sin}\theta .$$

Note that every 2-surface $`𝒮(t,r):(t=const,r=const)`$ is a homogeneous 2-sphere of area

$$A=4\pi S^2(t)f^2(r).$$

### 2.1 Radial Null Geodesics

The family of past-directed radial null geodesics in this space-time has tangent vector field

$$K^a=\frac{1}{S(t)}(1,\pm \frac{1}{S(t)},0,0)=\frac{dx^a}{dv}$$

 where the sign depends on whether the geodesics are ingoing or outgoing. This form gives

$$K^aK_a=g_{ab}K^aK^b=0,$$

$$K^au_a=\frac{1}{S(t)}=(1+z)$$

as required, and they are normal to the family of instantaneous homogeneous 2-spheres $`𝒮(t,r)`$ because they lie in the orthogonal 2-plane $`(\theta =const,\varphi =const)`$ to these 2-spheres, described by coordinates $`(t,r).`$ They diverge from the point of origin of coordinates $`r=0,`$ and refocus at the antipodal point $`r=\pi `$ when $`k=+1.`$

The divergence of this family of null geodesics is given by

$$K_{;a}^a=\frac{1}{\sqrt{g}}\frac{}{x^a}(\sqrt{g}K^a)=\frac{2}{S^2(t)}\left[\dot{S}(t)\pm \frac{f(r)/r}{f(r)}\right].$$
(4)
CTSs occur if the divergence is negative for both families of null geodesics for some value of $`r`$ and $`t`$, i.e. for both choices of sign in the last term on this 2-surface.We now look at a series of specific cases.

## 3 The Early and Late Universe

### 3.1 Standard radiation/matter dominated expansion

Here $`S(t)=at^n`$ with $`a>0`$ and $`n=2/3`$ in the matter era, $`n=1/2`$ in the radiation era, and with $`0<t<\mathrm{}`$. The divergence of the radial families of null geodesics is given by (4):

$$K_{;a}^a=\frac{2}{at^{2n}}\left[ant^{n1}\pm (\mathrm{cot}r,\frac{1}{r},\mathrm{coth}r)\right]$$

for $`k=+1,0.1`$ respectively. These expressions give the divergences of the normals of the 2-spheres $`(r,t)`$ constant. For each $`t>0`$ and each value of $`k`$ these will be CTSs, obtained by choosing $`r`$ large enough that the magnitude of the second term in the square brackets is less than that of the first term in these brackets, so that $`K_{;a}^a<0`$ for both signs (i.e. for both ingoing and outgoing null geodesics). Note that previous examinations of these surfaces have concentrated on showing that the past light cone of any observer will start reconverging and thus that closed trapped surfaces will occur in the past of the observer, associated with this refocussing of their past light cone. What is shown here is a bit different, namely that for any value of $`t,`$ if one goes to large enough values of $`r`$ there will be closed trapped surfaces surrounding the origin at that time. Hence there are closed trapped surfaces surrounding us even today (associated with the reconvergence of the past light cone of observers in our future).

Here the energy condition $`\mu +3p>0`$ is satisfied because we consider only ordinary matter. If this energy condition does not change in the past then a singularity is predicted in all cases via the Raychaudhuri equation for exact FL models , and via Theorems 2&3 above for universes that are perturbed FRW models at late times even if they are quie different at early times (small enough perturbations will necessarily preserve the inequalities $`K_{;a}^a<0`$ in the late-time era), and also by Theorem 1 in the case $`k1`$. Hence the existence of singularities is predicted and is stable to perturbations of these models at late times (which may correspond to very large changes in the models at early enough times) .

### 3.2 Inflation

In realistic universes the Hot Big Bang era may be preceded by an era of inflation. It will still be true that CTSs occur in the late universe (between decoupling and the present day, as well as in the Hot Big Bang era itself). However inflation violates the timelike convergence condition at early times and so can avoid the initial singularity that would otherwise be predicted because these CTSs exist at later times. Theorems 2 and 3 fall away because of the timelike energy condition violation. If $`k=0`$ or $`1`$, then a singularity will indeed occur, because of the Friedmann equation in the exact FL case , and in the perturbed FL case because of Theorem 1, relying only on the null energy condition, together with the existence of open space sections if they have their normal topology. Singularities are not inevitable when $`k=+1,`$ both because the Friedmann equation now allows a minimum and because there are then closed spacelike sections so none of the above theorems apply. Various kinds of non-singular model can then occur

## 4 The de Sitter Universe

We deal in turn with the three RW frames for de Sitter spacetime (see  for its global properties).

### 4.1 The k=+1 Frame

In the global $`k=+1`$ frame, the metric is (3) with $`S(t)=A\mathrm{cosh}Ht`$, $`f(r)=\mathrm{sin}r,`$ where $`A,`$ $`H>0`$ and $`\mathrm{}<t<\mathrm{}.`$ The cosmology is non-singular and geodesically complete. Note that $`0r\pi .`$ The antipodal point to the origin of coordinates is at $`r=\pi ;`$the equator for these coordinates is at $`r=\pi /2.`$ Then by (4)

$$K_{;a}^a=\frac{2}{A^2\mathrm{cosh}^2Ht}\left[AH\mathrm{sinh}Ht\pm \frac{\mathrm{cos}r}{\mathrm{sin}r}\right]$$

We get a past closed trapped surface if $`t>0`$ and $`AH\mathrm{sinh}Ht>|\frac{\mathrm{cos}r}{\mathrm{sin}r}|.`$ Now the latter term is zero at $`r=\pi /2`$ and diverges positively and negatively respectively at $`r=0,\pi .`$ Thus provided $`t>0,`$ we have a past closed trapped surface for all $`r`$ such that

$$AH\mathrm{sinh}Ht>|\mathrm{cot}r|$$

which defines a non-zero set of 2-spheres around the equator at $`r=\pi /2`$. As $`t0,`$ these shrink to just the equator; as $`t\mathrm{}`$ , they expand to a large part of the whole sphere. The origin is arbitrary, so every 2-sphere $`(t,r)`$ constant with $`t>0`$ and area greater than

$$A_{}=4\pi A^2\mathrm{cosh}^2Ht\mathrm{sin}^2(r_{}),AH\mathrm{sinh}Ht=\mathrm{cot}r_{}$$

in any $`k=+1`$ frame will be a past closed trapped surface. Its normals will self-intersect and have caustics where $`K_{;a}^a\mathrm{}`$ at both $`r0`$ (the origin, the ingoing family) and $`r\pi `$ (the antipodal point, the outgoing family). Since $`\mathrm{sin}^2r=\frac{1}{1+\mathrm{cot}^2r}`$ this gives

$$A_{}=4\pi \frac{A^2\mathrm{cosh}^2Ht}{1+A^2H^2\mathrm{sinh}^2Ht}.$$

As $`t0`$, $`S(t)A^2`$ and $`A_{}4\pi A^2`$ (corresponding to $`r\pi /2`$); as $`t\mathrm{}`$, $`S(t)\mathrm{}`$ and $`A_{}4\pi \frac{1}{H^2}`$, so there is a minimum radius that will give a CTS (but the size of the 3-spaces increases without limit, so this minimum radius will be an ever smaller fraction of the size of the space sections). For $`t<0`$ we find the corresponding family of future trapped surfaces. There are no trapped surface for $`t=0.`$

For these trapped surfaces (2) becomes

$`{\displaystyle \frac{d\theta }{dv}}`$ $`=`$ $`{\displaystyle \frac{1}{2}}\theta ^2{\displaystyle \frac{d\theta }{\theta ^2}}={\displaystyle \frac{1}{2}}dv{\displaystyle \frac{1}{\theta }}={\displaystyle \frac{1}{2}}(vv_0)`$
$``$ $`\theta ={\displaystyle \frac{2}{v}}\mathrm{on}\mathrm{choosing}v_0=0.`$
Thus the geodesics generating the pasts of the set of past trapped surfaces locally self-intersect, hence signalling an end to the boundary of the past of the trapped surfaces. But these intersections occur round the back (near the antipodal points); hence the past of these surfaces is not trapped by these null geodesics. Theorem 1 does not apply, even though the null energy condition is true, because the spatial sections are compact.

Note that these space sections are not unique : every 3-sphere passing through the de Sitter throat is equivalent to every other one. The same result as above must be found in every such frame. This at first seems to lead to an apparent contradiction: no consistency is found at the same points on the hyperboloid in different frames, because in the $`t=0`$ frames of different such choices there are no trapped surfaces but at the same points in the $`t0`$ frames there do exist closed trapped surfaces. However the correspondence between events on the hyperboloid and 2-spheres depends on the frame chosen. The frame corresponds point by point to events on the hyperboloid in each coordinate frame in a way such that the corresponding 2-spheres (spherically surrounding the origin of coordinates) depends on the coordinate frame chosen. Given a specific choice of frame, however, a unique such correspondence of points and 2-spheres exists . Those found in the $`t=0`$surfaces in a particular frame may be represented as marginally trapped 2-spheres in that frame, but are fully trapped in other frames. Thus one must choose a specific frame and work it all out in that frame; the results in all other frames will then follow by boosting, rotating and translating that frame.

### 4.2 Perturbed de Sitter Universes

The further basic point is that existence of CTSs do not imply singularities in perturbed de Sitter universes either, when the null generic condition holds, even when $`\mu +p0`$ (so that self-intersections occur because of the generic conditionds as outlined above) but provided still $`\mu +3p<0`$. This is because they don’t imply them in the de Sitter case, where the past of each CTS is also compact ,but that does not imply a singularity because of the closed space sections. If this were not true the de Sitter universe would be unstable - but it is well known to be stable.

### 4.3 The k=0 frame

In the $`k=0`$ frame, which covers half the spacetime (and so is not geodesically complete), we have the metric (3) with $`S(t)=A\mathrm{exp}Ht`$ , $`f(r)=r,`$ for $`\mathrm{}<t<\mathrm{}`$ and $`A,H>0.`$ In this case (4 ) shows

$$K_{;a}^a=\frac{2}{A^2\mathrm{exp}2Ht}\left[AH\mathrm{exp}Ht\pm \frac{1}{r}\right]$$

The second term dies away to zero, so there will be closed trapped surfaces for

$$AH\mathrm{exp}Ht>\frac{1}{r}$$

which will always be true for large enough $`r`$ for any $`t.`$ These will be part of the same set of 2-spheres as characterised above, but expressed in different coordinates. In this case they do correspond to geodesic incompleteness, because this coordinate frame does not cover the whole hyperboloid and Theorem 1 above applies in this frame. New information could come in from the other half of the hyperboloid if the solution is extended further (which even though it lies beyond the infinite redshift surface need not correspond to infinite redshift of matter beyond that surface; that would depend on how matter is moving in this further part of the hyperboloid). The spacetime is singular but extendible.

### 4.4 The k=-1 frame

In the $`k=1`$ frame, which covers less than half the spacetime (and so is not geodesically complete), we have the metric (3) with $`S(t)=A\mathrm{sinh}Ht`$, $`f(r)=\mathrm{sinh}r,`$ for $`0<t<\mathrm{}`$ and $`A,H>0.`$ In this case (4) shows that

$$K_{;a}^a=\frac{2}{A^2\mathrm{sinh}^2Ht}\left[AH\mathrm{cosh}Ht\pm \frac{\mathrm{cosh}r}{\mathrm{sinh}r}\right]$$

The latter term is always of magnitude $`>1`$, diverging as $`r0`$ and $`1`$ as $`r\mathrm{}.`$ Thus there will be closed trapped surface for values of time $`t`$ such that

$$AH\mathrm{cosh}Ht>1,t>0t>(1/H)\mathrm{arg}\mathrm{cosh}(1/AH)$$

which will exist for all $`A,H>0.`$ For those values of $`t`$, closed trapped surfaces exist for all values of $`r`$ such that

$$AH\mathrm{cosh}Ht>\mathrm{coth}r$$

which will then exist for large enough $`r.`$ These will again be part of the same set of 2-spheres as characterised above, but expressed in different coordinates. Here again Theorem 1 predicts an in initial singularity and again they are geodesically incomplete but extendible.

#### 4.4.1 Conclusion: de Sitter spacetime

In all three cases, we show that past closed trapped surfaces exist in the de Sitter universe (of course these are just the same set of 2-surfaces found in different coordinates) and lead to the boundaries of the pasts of those 2-surfaces being compact. This does not however lead to spacetime singularities in the first case (de Sitter spacetime is geodesically complete in the $`k=+1`$ frame), because the pasts of the trapped 2-spheres are not trapped by these null boundaries, rather they can escape freely to earlier times because the 3-spaces are compact (the null cone intersections take place at the antipodal point on the other side of the 3-spaces with 3-sphere topology, allowing any interior matter to escape on this side, near the origin). However in the $`k=0`$ and $`k=1`$ frames the cosmology is singular because the worldlines do not cover the whole spacetime; indeed they self-intersect at $`t=0`$ in the $`k=1`$ case.

## 5 Emergent Universes

These are non-singular models with $`k=+1`$ that start off asymptotically as Einstein Static universes, and then evolve to de Sitter universes (and at even later times to a shandard hot big bang) . Here $`k=+1`$ and in the simplest case $`S(t)=A+B\mathrm{exp}Ht`$, with $`\mathrm{}<t<\mathrm{},`$ and $`A,B>0`$ . Then by (4)

$$K_{;a}^a=\frac{2}{(A+B\mathrm{exp}Ht)^2}\left[BH\mathrm{exp}Ht\pm \frac{\mathrm{cos}r}{\mathrm{sin}r}\right]$$

We get a past closed trapped surface if $`BH\mathrm{exp}Ht>|\frac{\mathrm{cos}r}{\mathrm{sin}r}|.`$ Now the latter is zero at $`r=\pi /2`$ and diverges to -$`\mathrm{},\mathrm{}`$ at $`r=0,\pi `$respectively.Thus provided $`t>0,`$ we have a closed trapped surface for all $`r`$ such that

$$BH\mathrm{exp}Ht>|\mathrm{cot}r|$$

which defines a non-zero set of 2-spheres around the equator, depending on the coordinate time $`t`$. As $`t\mathrm{},`$ the static limit, these shrink to just the equator; as $`t\mathrm{}`$, they expand to a large part of the whole sphere (as in the de Sitter case). However the spacetime is geodesically complete and non-singular. The past generators of the 2-spheres intersect , and this does not imply the exstence of singularities; none of the singularity theorems apply. The de Sitter universe is the special case when $`A=0,`$ obtained as a completely smooth limit.

### 5.1 Perturbed Emergent universes

The existence of CTSs does not imply singularities in perturbed emergent universe eithers. Perturbations that ensure the genericity condition on all null geodesics again do not imply a singularity, as in the case of the de Sitter universe.

## 6 The Einstein Static Universe

Here we have a RW metric (3) with $`S(t)=S_0>0,`$ $`k=+1,`$ $`\mathrm{}<t<\mathrm{}`$, and so

$$K_{;a}^a=\frac{2}{S^2(t)}\left[\dot{S}(t)\pm \frac{f/r}{f(r)}\right]=\frac{2}{S_0^2}\left[0\pm \frac{\mathrm{cos}r}{\mathrm{sin}r}\right].$$

In this case there are no closed trapped surfaces; however marginally trapped surfaces occur on the equator, the null geodesics generating their pasts intersecting at the antipodal pointRoughly: the closure of the space sections does not allow existence of 2-spheres that are large enough to be trapped. . This lack of closed trapped surfaces is connected to the high degree of stability of the E-S universe (they are stable to all inhomogeneus perturbations in the radiation case ). Again it is thus true that perturbing the universe leads to the generic condition for null geodesics but not to singularities.

## 7 Conclusion

Closed trapped surfaces occur in most Friedmann models, including the de Sitter universe. They necessarily lead to a singularity only if $`\rho +3p>0`$ . When $`\rho +3p<0`$ and $`k=+1`$, singularity avoidance is possible. The null energy condition $`\rho +p>0`$ does not necessarily lead to a singularity, despite existence of these closed trapped surfaces and hence compact boundaries of the past of these 2-spheres, when the spatial sections are compact.


format: text

text:
# On the foliation of space-time by constant mean curvature hypersurfaces

## 0. Introduction

In  it is proved that a globally hyperbolic Lorentzian manifold $`N`$ with a compact Cauchy hypersurface can be foliated by constant mean curvature hypersurfaces if the big bang and  big crunch hypotheses are valid, and if the time-like convergence condition holds, i.e. if

(0.1)
$$\overline{R}_{\alpha \beta }\eta ^\alpha \eta ^\beta 0$$

for all time-like vectorfields $`(\eta ^\alpha )`$.

If we assume for simplicity that the mean curvature of the barriers provided by the big bang resp. big crunch hypotheses tend to $`\mathrm{}`$ resp. $`+\mathrm{}`$, the foliation can be described as follows: Indicate by $`M_\tau `$ a closed hypersurface of mean curvature $`\tau `$, then, the foliation consists of the uniquely determined $`M_\tau `$, $`0\tau `$, and of the set $`𝒞_0`$ of maximal slices. If there is more than one maximal slice, then $`𝒞_0`$ comprises a whole continuum of maximal slices, which are all totally geodesic and the ambient metric is static in $`𝒞_0`$.

The mean curvature of the slices of the foliation can be looked at as a function $`\tau `$ on $`N`$, which is continuous as one easily checks. However, in order to use $`\tau `$ as a new time function, $`\tau `$ has to be smooth with non-vanishing gradient.

We prove that this is indeed the case in $`\{\tau 0\}`$ and also globally, if there is a maximal slice that is not totally geodesic or if the strict inequality is valid in (0.1). Evidently, if we have two (and more) maximal slices, then, $`D\tau `$ vanishes in the interior of $`𝒞_0`$, and we shall also give an example of a foliation with exactly one totally geodesic slice, where $`D\tau `$ vanishes on that slice, i.e. the assumptions guaranteeing $`\tau `$ to be smooth and $`D\tau 0`$ are also sharp.

The paper is organized as follows: In Section 1 we introduce the notations and common definitions we rely on.

The main theorem is proved in Section 2, while the counterexample is given in Section 3.

## 1. Notations and definitions

The main objective of this section is to state the equations of Gauß, Codazzi, and Weingarten for space-like hypersurfaces $`M`$ in a (n+1)-dimensional Lorentzian manifold $`N`$. Geometric quantities in $`N`$ will be denoted by $`(\overline{g}_{\alpha \beta }),(\overline{R}_{\alpha \beta \gamma \delta })`$, etc., and those in $`M`$ by $`(g_{ij}),(R_{ijk)}`$, etc. Greek indices range from $`0`$ to $`n`$ and Latin from $`1`$ to $`n`$; the summation convention is always used. Generic coordinate systems in $`N`$ resp. $`M`$ will be denoted by $`(x^\alpha )`$ resp. $`(\xi ^i)`$. Covariant differentiation will simply be indicated by indices, only in case of possible ambiguity they will be preceded by a semicolon, i.e. for a function $`u`$ in $`N`$, $`(u_\alpha )`$ will be the gradient and $`(u_{\alpha \beta })`$ the Hessian, but e.g., the covariant derivative of the curvature tensor will be abbreviated by $`\overline{R}_{\alpha \beta \gamma \delta ;ϵ}`$. We also point out that

(1.1)
$$\overline{R}_{\alpha \beta \gamma \delta ;i}=\overline{R}_{\alpha \beta \gamma \delta ;ϵ}x_i^ϵ$$

with obvious generalizations to other quantities.

Let $`M`$ be a space-like hypersurface, i.e. the induced metric is Riemannian, with a differentiable normal $`\nu `$ which is time-like.

In local coordinates, $`(x^\alpha )`$ and $`(\xi ^i)`$, the geometric quantities of the space-like hypersurface $`M`$ are connected through the following equations

(1.2)
$$x_{ij}^\alpha =h_{ij}\nu ^\alpha $$

the so-called Gauß formula. Here, and also in the sequel, a covariant derivative is always a full tensor, i.e.

(1.3)
$$x_{ij}^\alpha =x_{,ij}^\alpha \mathrm{\Gamma }_{ij}^kx_k^\alpha +\overline{\mathrm{\Gamma }}_{\beta \gamma }^\alpha x_i^\beta x_j^\gamma .$$

The comma indicates ordinary partial derivatives.

In this implicit definition the second fundamental form $`(h_{ij})`$ is taken with respect to $`\nu `$.

The second equation is the Weingarten equation

(1.4)
$$\nu _i^\alpha =h_i^kx_k^\alpha ,$$

where we remember that $`\nu _i^\alpha `$ is a full tensor.

Finally, we have the Codazzi equation

(1.5)
$$h_{ij;k}h_{ik;j}=\overline{R}_{\alpha \beta \gamma \delta }\nu ^\alpha x_i^\beta x_j^\gamma x_k^\delta $$

and the Gauß equation

(1.6)
$$R_{ijkl}=\{h_{ik}h_{jl}h_{il}h_{jk}\}+\overline{R}_{\alpha \beta \gamma \delta }x_i^\alpha x_j^\beta x_k^\gamma x_l^\delta .$$

Now, let us assume that $`N`$ is a globally hyperbolic Lorentzian manifold with a compact Cauchy surface. $`N`$ is then a topological product $`\times 𝒮_0`$, where $`𝒮_0`$ is a compact Riemannian manifold, and there exists a Gaussian coordinate system $`(x^\alpha )`$, such that the metric in $`N`$ has the form

(1.7)
$$d\overline{s}_N^2=e^{2\psi }\{dx_{}^{0}{}_{}{}^{2}+\sigma _{ij}(x^0,x)dx^idx^j\},$$

where $`\sigma _{ij}`$ is a Riemannian metric, $`\psi `$ a function on $`N`$, and $`x`$ an abbreviation for the space-like components $`(x^i)`$, see , \[7, p. 212\], \[6, p. 252\], and \[3, Section 6\]. We also assume that the coordinate system is future oriented, i.e. the time coordinate $`x^0`$ increases on future directed curves. Hence, the contravariant time-like vector $`(\xi ^\alpha )=(1,0,\mathrm{},0)`$ is future directed as is its covariant version $`(\xi _\alpha )=e^{2\psi }(1,0,\mathrm{},0)`$.

Let $`M=\mathrm{graph}u_{|_{𝒮_0}}`$ be a space-like hypersurface

(1.8)
$$M=\{(x^0,x):x^0=u(x),x𝒮_0\},$$

then the induced metric has the form

(1.9)
$$g_{ij}=e^{2\psi }\{u_iu_j+\sigma _{ij}\}$$

where $`\sigma _{ij}`$ is evaluated at $`(u,x)`$, and its inverse $`(g^{ij})=(g_{ij})^1`$ can be expressed as

(1.10)
$$g^{ij}=e^{2\psi }\{\sigma ^{ij}+\frac{u^i}{v}\frac{u^j}{v}\},$$

where $`(\sigma ^{ij})=(\sigma _{ij})^1`$ and

(1.11) $`u^i`$ $`=\sigma ^{ij}u_j`$
$`v^2`$ $`=1\sigma ^{ij}u_iu_j1|Du|^2.`$
Hence, $`\mathrm{graph}u`$ is space-like if and only if $`|Du|<1`$.

The covariant form of a normal vector of a graph looks like

(1.12)
$$(\nu _\alpha )=\pm v^1e^\psi (1,u_i).$$

and the contravariant version is

(1.13)
$$(\nu ^\alpha )=v^1e^\psi (1,u^i).$$

Thus, we have

###### Remark 1.1.

Let $`M`$ be space-like graph in a future oriented coordinate system. Then, the contravariant future directed normal vector has the form

(1.14)
$$(\nu ^\alpha )=v^1e^\psi (1,u^i)$$

and the past directed

(1.15)
$$(\nu ^\alpha )=v^1e^\psi (1,u^i).$$

In the Gauß formula (1.2) we are free to choose the future or past directed normal, but we stipulate that we always use the past directed normal for reasons that we have explained in \[4, Section 2\].

Look at the component $`\alpha =0`$ in (1.2) and obtain in view of (1.15)

(1.16)
$$e^\psi v^1h_{ij}=u_{ij}\overline{\mathrm{\Gamma }}_{00}^0u_iu_j\overline{\mathrm{\Gamma }}_{0j}^0u_i\overline{\mathrm{\Gamma }}_{0i}^0u_j\overline{\mathrm{\Gamma }}_{ij}^0.$$

Here, the covariant derivatives are taken with respect to the induced metric of $`M`$, and

(1.17)
$$\overline{\mathrm{\Gamma }}_{ij}^0=e^\psi \overline{h}_{ij},$$

where $`(\overline{h}_{ij})`$ is the second fundamental form of the hypersurfaces $`\{x^0=\text{const}\}`$.

An easy calculation shows

(1.18)
$$\overline{h}_{ij}e^\psi =\frac{1}{2}\dot{\sigma }_{ij}\dot{\psi }\sigma _{ij},$$

where the dot indicates differentiation with respect to $`x^0`$.

## 2. $`\tau `$ is a time function

Let $`M_0=M_{\tau _0}`$ be a hypersurface of constant mean curvature $`\tau _0`$, and let $`(x^\alpha )`$ be a future oriented, normal Gaussian coordinate system relative to $`M_0`$, i.e.

(2.1)
$$d\overline{s}_N^2=dx_{}^{0}{}_{}{}^{2}+\sigma _{ij}(x^0,x)dx^idx^j$$

in a tubular neighbourhood $`𝒰`$ of $`M_0`$; here, $`x`$ is an abbreviation for the spatial components $`(x^i)`$.

Since we have a continuous foliation of $`N`$, the constant mean curvature slices contained in $`𝒰`$ can be written as graphs over $`M_0`$

(2.2)
$$M_\tau =\{(x^0,x^i):x^0=u(\tau ,x^i)\}$$

with a continuous function $`uC^0(I\times M_0)`$, $`I`$ an open interval containing $`\tau _0`$.

Let us assume for the moment that $`u`$ is smooth—we shall see later that this is indeed the case, at least for $`\tau _00`$—, and let us define the transformation

(2.3)
$$\mathrm{\Phi }(\tau ,x^i)=(u(\tau ,x^i),x^i).$$

Then, we have

(2.4)
$$detD\mathrm{\Phi }=\frac{u}{\tau }\dot{u}.$$

In view of the monotonicity of the foliation, cf. \[3, Lemma 7.2\], we know that $`\dot{u}0`$. Thus, if we could show that $`u`$ is smooth and $`\dot{u}`$ strictly positive, we would obtain that $`\mathrm{\Phi }`$ is a diffeomorphism, and hence, that $`\tau `$ is smooth with non-vanishing gradient.

Let us first show that $`u`$ is smooth.

###### Lemma 2.1.

Let $`\tau _0`$ be such that $`M_{\tau _0}`$ is not totally geodesic or assume that the strict inequality is valid in (0.1). Consider a tubular neighbourhood $`𝒰`$ of $`M_0=M_{\tau _0}`$, and choose normal Gaussian coordinates as above, such that the $`M_\tau 𝒰`$ are graphs

(2.5)
$$M_\tau =\mathrm{graph}u(\tau )_{|_{M_0}}.$$

Then, $`u`$ is smooth in $`I\times M_0`$, where $`I=(\tau _0ϵ,\tau _0+ϵ)`$ for small $`ϵ=ϵ(\tau _0)>0`$.

###### Proof.

Let $`u_0=u(\tau _0)`$, which is incidentally identical zero due to our choice of coordinates, and define the operator $`G`$

(2.6)
$$G(\tau ,\phi )=H(\phi )\tau ,$$

where $`H(\phi )`$ is an abbreviation for the mean curvature of $`\mathrm{graph}\phi _{|_{M_0}}`$.

It is well-known, cf. \[2, Section 5\], that $`D_2G(\tau _0,u_0)\phi `$ equals

(2.7)
$$\mathrm{\Delta }\phi +\phi \{A^2+\overline{R}_{\alpha \beta }\nu ^\alpha \nu ^\beta \},$$

where $`A^2=h_{ij}h^{ij}`$, and the Laplacian, the second fundamental form and $`\nu `$ are evaluated with respect to $`M_0`$.

Hence, $`D_2G`$ is invertible at $`(\tau _0,u_0)`$ in view of our assumptions, and thus, the implicit function theorem is applicable to yield the existence of a smooth function $`u=u(\tau ,x)`$ satisfying

(2.8)
$$G(\tau ,u(\tau ,x))=0$$

for $`|\tau \tau _0|<ϵ`$, $`ϵ`$ small; of course, the new smooth function coincides with the original continuous function. ∎

It remains to demonstrate that $`\dot{u}`$ is strictly positive.

We use the same setting as above, i.e. we work in a tubular neighbourhood $`𝒰`$ of a slice $`M_{\tau _0}`$ with corresponding normal Gaussian coordinates, and look at the $`M_\tau 𝒰`$, which are then graphs over $`M_{\tau _0}`$.

###### Lemma 2.2.

Let $`M_\tau =\mathrm{graph}u(\tau )`$ and $`M_{\overline{\tau }}=\mathrm{graph}u(\overline{\tau })`$, then, there exists a constant $`c=c(𝒰)`$ such that

(2.9)
$$|\tau \overline{\tau }|c\underset{M_{\tau _0}}{inf}|u\overline{u}|,$$

where $`\overline{u}=u(\overline{\tau })`$.

###### Proof.

Assume for simplicity $`\overline{\tau }<\tau `$. We want to employ the relation (1.16)—now with $`\psi 0`$—, but the covariant derivatives of $`u`$ should be expressed with respect to the metric $`(\sigma _{ij}(u,x))`$; recall that

(2.10)
$$g_{ij}=u_iu_j+\sigma (u,x).$$

An easy calculation reveals that

(2.11)
$$u_{ij}=v^2u_{;ij},$$

where the semicolon indicates the covariant derivatives with respect to the metric $`(\sigma _{ij}(u,x))`$. Let $`\stackrel{~}{\mathrm{\Gamma }}_{ij}^k`$ be the corresponding Christoffel symbols and indicate by a comma ordinary partial derivatives, then, we obtain

(2.12)
$$h_{ij}=v^1u_{,ij}+v^1\stackrel{~}{\mathrm{\Gamma }}_{ij}^ku_kv\overline{\mathrm{\Gamma }}_{ij}^0$$

with an analogous relation valid for $`\overline{u}`$; note that the other Christoffel symbols in (1.16) vanish, since $`\psi 0`$.

Let $`x_0M_{\tau _0}`$ be such that

(2.13)
$$(u\overline{u})(x_0)=\underset{M_{\tau _0}}{inf}(u\overline{u}).$$

We now observe that in $`x_0`$

(2.14)
$$Du=D\overline{u}\text{and}D^2uD^2\overline{u},$$

and deduce from (2.12), and the corresponding relation for $`M_{\overline{\tau }}`$,

(2.15)
$$\tau \overline{\tau }f(x_0,u(x_0))f(x_0,\overline{u}(x_0))$$

with a smooth function $`f`$, i.e.

(2.16)
$$\tau \overline{\tau }c[u(x_0)\overline{u}(x_0)].$$

###### Remark 2.3.

An estimate of the form (2.9) is valid in any Gaussian coordinate system in which the slices of the foliation are represented as graphs. The constant $`c`$ depends on the compact set containing the slices under consideration.

Thus, we have proved

###### Theorem 2.4.

The function $`\tau `$ is smooth with non-vanishing gradient in $`\{\tau 0\}`$. If there is a maximal slice that is not totally geodesic or if the strict inequality is valid in (0.1), then $`\tau `$ is a globally defined time function.

We note that by assumption there are maximal slices in the foliation.

## 3. A counterexample

Let $`S^n`$ be the standard unit sphere with metric $`\sigma _{ij}`$, and for $`tI=(ϵ,ϵ)`$ let

(3.1)
$$f=f(t)=_0^t\frac{s^3}{ϵ^2s^2}𝑑s.$$

Then, we define $`N=I\times S^n`$ with metric

(3.2)
$$d\overline{s}_N^2=dt^2+e^{2f}\sigma _{ij}dx^idx^j,$$

where the coordinates $`(t,x^i)`$ are supposed to be future oriented.

The level hypersurfaces $`\{t=\text{const}\}`$ are all totally umbilical, and their mean curvature is equal to

(3.3)
$$\tau =H=n\dot{f}=n\frac{t^3}{ϵ^2t^2}\text{.}$$

Thus, we see that the big bang and big crunch hypotheses are satisfied and that $`D\tau `$ vanishes on the unique totally geodesic maximal slice.

It remains to verify the time-like convergence condition.

We first note, that for $`ϵ1`$

(3.4)
$$\ddot{f}+\dot{f}^2=\frac{3t^2}{ϵ^2t^2}\frac{2t^4}{(ϵ^2t^2)^2}+\frac{t^6}{(ϵ^2t^2)^2}0,$$

which will be the core inequality to estimate the Ricci tensor.

Secondly, let us write the metric as a conformal metric

(3.5)
$$\overline{g}_{\alpha \beta }=e^{2\psi }g_{\alpha \beta }$$

in coordinates $`(x^0,x^i)`$, where

(3.6) $`g_{\alpha \beta }dx^\alpha dx^\beta `$ $`=dx_{}^{0}{}_{}{}^{2}+\sigma _{ij}dx^idx^j,`$
(3.7) $`{\displaystyle \frac{dx^0}{dt}}`$ $`=e^{f(t)},`$
and $`\psi (x^0)=f(t)`$.

The metric $`(g_{\alpha \beta })`$ is a product metric and, therefore, the only non-zero components of its Ricci tensor $`(R_{\alpha \beta })`$ are of the form $`R_{ij}`$ and coincide with the components of the Ricci tensor of $`S^n`$.

The Ricci tensor $`(\overline{R}_{\alpha \beta })`$ of the metric $`(\overline{g}_{\alpha \beta })`$ is connected to $`(R_{\alpha \beta })`$ through the relation

(3.8)
$$\begin{array}{cc}\hfill \overline{R}_{\alpha \beta }& =R_{\alpha \beta }(n1)[\psi _{\alpha \beta }\psi _\alpha \psi _\beta ]\hfill \\ & g_{\alpha \beta }[\mathrm{\Delta }\psi +(n1)|D\psi |^2],\hfill \end{array}$$

where covariant derivatives and the norm are calculated with respect to $`(g_{\alpha \beta })`$; note that $`|D\psi |^2`$ can be negative.

Hence, we conclude

(3.9) $`\overline{R}_{00}`$ $`=n\ddot{\psi }=ne^{2f}[\ddot{f}+\dot{f}^2],`$
$`\overline{R}_{0i}`$ $`=0,`$
and

(3.10) $`\overline{R}_{ij}`$ $`=R_{ij}+\sigma _{ij}[\ddot{\psi }+(n1)\dot{\psi }^2],`$
$`=R_{ij}+\sigma _{ij}[\ddot{f}+n\dot{f}^2]e^{2f}.`$
Let $`(\eta ^\alpha )`$ be a time-like vector with component $`\eta ^0=1`$, so that

(3.11)
$$\sigma _{ij}\eta ^i\eta ^j<1.$$

Then, taking into account that $`\ddot{f}0`$, we derive

(3.12)
$$\begin{array}{cc}\hfill \overline{R}_{\alpha \beta }\eta ^\alpha \eta ^\beta & =ne^{2f}[\ddot{f}+\dot{f}^2]+R_{ij}\eta ^i\eta ^j\hfill \\ & +\sigma _{ij}\eta ^i\eta ^j[\ddot{f}+n\dot{f}^2]e^{2f}\hfill \\ & ne^{2f}[\ddot{f}+\dot{f}^2][1+\sigma _{ij}\eta ^i\eta ^j],\hfill \end{array}$$

and we conclude that the time-like convergence condition is satisfied in view of (3.4).


text:
# Numerical 3D constraints on convective eddy time-correlations : consequences for stochastic excitation of solar p modes.

## 1 Introduction

Solar oscillations are believed to be stochastically excited by turbulent convection in the outer part of the Sun. The excitation is caused by turbulent convective motions which generate acoustic energy which in turn is injected into the p modes.

Models of stochastic excitation of stellar p modes have been proposed by several authors, (e.g. Goldreich & Keeley, 1977; Osaki, 1990; Balmforth, 1992; Goldreich et al., 1994). These models use simplified models to describe the dynamics of the turbulent medium. For instance these approaches (Goldreich & Keeley, 1977; Balmforth, 1992) assume a gaussian function for representing $`\chi _k`$, the frequency component of the correlation product of the turbulent velocity field. As pointed out by Samadi (2001), the way the component $`\chi _k`$ is modeled plays a crucial role in controlling the extent of the excitation region of a given mode and hence the total amount of acoustic energy injected into the mode. In the following, $`\chi _k`$ will also be referred to as *the dymamic model of turbulence* and *dynamic* will refer to time-dependence or frequency-dependence.

Direct computations of the rate at which the solar p modes are excited have been performed by Stein & Nordlund (2001) using 3D simulations of the upper part of the solar convective zone. They found good agreement between their numerical results and the solar seismic observations. This direct but time consuming approach did not address the role of the dynamic properties of the turbulent medium on the excitation mechanism.

In contrast semi-analytical formulations for $`P(\nu )`$ offer the advantage of testing separately several properties entering the excitation mechanism. Here we consider the formulation by Samadi & Goupil (2001, hereafter Paper I, see also () for a summary) which includes a detailed treatment of the *time averaged* and *dynamic* properties of the turbulent convective medium.

The impact of the averaged properties have been investigated by Samadi et al. (2003). The authors used a 3D simulation of the upper part of the solar convective zone to constrain the averaged properties of the turbulent convective medium. The computed rates $`P`$ at which the solar p modes are excited were found to be larger than those computed with a 1D mixing-length solar model but still *underestimate* the solar seismic data by a factor $`2.5`$. It was also found that the Reynolds tensor contributes about $`20\%`$ of the total acoustic energy injected into the solar p modes, in contrast with direct 3D estimations (Stein & Nordlund, 2001). These discrepancies were attributed to the assumed gaussian function for the dynamic model of turbulence.

In the present paper we therefore derive an empirical dynamic model of turbulence obtained from a 3D simulation of the upper part of the solar convective zone and then study the consequences of using this model on the computed excitation rates $`P`$. We compare our computation with solar seismic data and finally obtain on an improved model of stochastic excitation.

The paper is organised as follows: The basic theoretical background and notations are recalled in Sect. 2. In Sect. 3, a 3D simulation of the upper part of the solar convective zone is used to characterise $`\chi _k`$ in the domain where stochastic excitation takes place. The inferred $`\nu `$-dependency of $`\chi _k`$ is compared with the gaussian function and fitted with different non-gaussian functions. These functions are then used in Sect. 4 to compute the excitation rate $`P`$ for radial p modes. The results are compared with solar seismic observations as provided by Chaplin et al. (1998) and with computations in which the gaussian function is assumed. Sect. 5 is dedicated to discussions and conclusions.

## 2 Stochastic excitation

### 2.1 The model of stochastic excitation

We consider the model of stochastic excitation as described in Paper I and assume here - as in Samadi et al. (2003) \- that injection of acoustic energy into the modes is isotropic and consider only radial p modes. Accordingly, the rate at which a given mode with frequency $`\omega _0`$ is excited can be written as :

$`P(\omega _0)`$ $`=`$ $`{\displaystyle \frac{\pi ^3}{2I}}{\displaystyle _0^M}\mathrm{d}m{\displaystyle \frac{\mathrm{\Phi }}{3}}\rho _0w^4\{{\displaystyle \frac{16}{15}}{\displaystyle \frac{\mathrm{\Phi }}{3}}\left({\displaystyle \frac{d\xi _\mathrm{r}}{dr}}\right)^2S_R`$ (1)
$`+{\displaystyle \frac{4}{3}}\left({\displaystyle \frac{\alpha _s\stackrel{~}{s}}{\rho _0w}}\right)^2{\displaystyle \frac{g_\mathrm{r}}{\omega _0^2}}S_S\}.`$
In Eq. (1), $`\rho _0`$ is the mean density, $`\xi _\mathrm{r}`$ is the radial component of the fluid displacement adiabatic eigenfunction $`\xi `$, $`I`$ is the mode inertia (Eq. (19)), $`\alpha _s=\left(p/s\right)_\rho `$ where $`p`$ denotes the pressure and $`s`$ the entropy, $`\stackrel{~}{s}`$ is the rms value of the entropy fluctuations which are assumed to arise solely from turbulence, $`g_r(\xi _\mathrm{r},r)`$ is a function involving the first and the second derivatives of $`\xi _r`$ with respect to $`r`$, $`\mathrm{\Phi }`$ is a mean anisotropy factor defined by Gough (1977) as

$$\mathrm{\Phi }(r)\frac{\overline{<𝐮^2><𝐮>^2}}{w^2(r)}$$
(2)
where $`𝐮`$ is the velocity field, $`<.>`$ denotes horizontal average, $`\overline{()}`$ denotes time average, and $`w(r)`$ is the mean vertical velocity ($`w^2\overline{<u_z^2><u_z>^2}`$). Expressions for $`g_r(\xi _\mathrm{r},r)`$ are given in Samadi et al. (2003).

The driving sources $`S_R(r,\omega _0)`$ and $`S_S(r,\omega _0)`$ arise from the Reynolds stress and the entropy fluctuations respectively :

$`S_R(r,\omega _0)`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}{\displaystyle \frac{\mathrm{d}k}{k^2}}{\displaystyle \frac{E(k,r)}{u_0^2}}{\displaystyle \frac{E(k,r)}{u_0^2}}\chi _k(\omega _0,r)`$ (3)
$`S_S(r,\omega _0)`$ $`=`$ $`{\displaystyle _0^{\mathrm{}}}{\displaystyle \frac{\mathrm{d}k}{k^2}}{\displaystyle \frac{E(k,r)}{u_0^2}}{\displaystyle \frac{E_s(k,r)}{\stackrel{~}{s}^2}}`$ (4)
$`\times {\displaystyle _{\mathrm{}}^+\mathrm{}}\mathrm{d}\omega \chi _k(\omega _0+\omega ,r)\chi _k(\omega ,r)`$
where $`u_0(r)\sqrt{\mathrm{\Phi }/3}w`$ is introduced for convenience, $`E(k,r)`$ is the time averaged turbulent kinetic energy spectrum, $`E_s(k,r)`$ is the time averaged turbulent spectrum associated with the entropy fluctuations and $`\chi _k(\omega ,r)`$ is the frequency-dependent part of the correlation product of the turbulent velocity field (see Section 2.2). In order to simplify the notation, we drop the explicit $`r`$ dependence of the quantities in Eqs. (1-4).

### 2.2 The dynamic model of turbulence

The dynamic model of turbulence is represented by $`\chi _k(\omega )`$. In order to give a precise meaning to $`\chi _k(\omega )`$, we recall first some theoretical relations. Excitation by Reynolds stresses involves $`\varphi _{i,j}(k,\omega )`$, the Fourier transform of the second-order velocity correlations; here the indices $`i`$ and $`j`$ refer to any of the 3 directions of the velocity field. For incompressible, homogeneous and isotropic turbulence, $`\varphi _{ij}(k,\omega )`$ has the form (Batchelor, 1970) :

$`\varphi _{ij}(k,\omega )`$ $`=`$ $`{\displaystyle \frac{E(k,\omega )}{4\pi k^2}}\left(\delta _{ij}{\displaystyle \frac{k_ik_j}{k^2}}\right)`$ (5)
where $`E(k,\omega )`$ is the turbulent kinetic energy spectrum as a function of $`k`$ and $`\omega `$ and $`\delta _{ij}`$ is the Kronecker tensor. Following Stein (1967), $`E(k,\omega )`$ is decomposed as

$`E(k,\omega )=E(k)\chi _k(\omega )`$ (6)
where $`\chi _k(\omega )`$ satisfies the normalisation condition (Tennekes & Lumley, 1982, Chap 8.1) :

$$_{\mathrm{}}^+\mathrm{}𝑑\omega \chi _k(\omega )=1.$$
(7)
According to the decomposition of Eq. (6) , $`\chi _k(\omega )`$ is - at fixed $`k`$ \- the frequency component of $`E(k,\omega )`$.

According to Eq. (5) and (6), $`\chi _k(\omega )`$ then represents the frequency dependence of $`\varphi _{i,j}(k,\omega )`$. In other words, $`\chi _k(\omega )`$ measures - in the frequency and $`k`$ wavenumber domains - the evolution of the velocity correlation between two distant points of the turbulent medium.

The same decomposition of Eq. (7) is assumed for $`E_s(k,\omega )`$. This leads to introducing $`\chi _k^s`$, the frequency-dependent part of the correlation product of the entropy fluctuation. For simplifying the computation of $`P`$, as in Paper I, we assume $`\chi _k^s=\chi _k`$. We have checked that $`\chi _k^s`$ and $`\chi _k`$ have almost the same behaviour in the region where excitation by entropy fluctuations is significant.

In the present work, we consider only the excitation of radial p modes. Let $`E_z(k,\omega )`$ be the vertical component of the kinetic energy spectrum. We consider that $`E_z(k,\omega )`$ can be decomposed as $`E(k,\omega )`$ (Eq. (7)). For isotropic turbulence we then have $`E(k,\omega )=3E_z(k,\omega )`$, $`E(k)=3E_z(k)`$ and $`\chi _k^z=\chi _k`$, which is equivalent to stating that the averaged and dynamic properties of the velocity field are the same in all 3 directions.

The anisotropy factor $`\mathrm{\Phi }`$ \- introduced in the expression for $`P`$, Eq. (1) - partially takes into account the spatial and temporal anisotropy of the turbulence ($`\mathrm{\Phi }=3`$ corresponds to a isotropic turbulence). It has been found in Samadi et al. (2003) that $`\mathrm{\Phi }2`$ within the region where most of the excitation occurs. This shows that the time and space averaged properties of the medium are indeed anisotropic. One can therefore expect that the dynamic properties of turbulence differ between the horizontal and vertical directions. As the excitation of radial p modes is predominantly governed by turbulent elements moving in the vertical direction, an open question is whether one should consider $`\chi _k^z`$ rather than $`\chi _k`$ in Eq. (3) and Eq. (4) when taking into account the dynamic properties of the turbulence. In the present work, we therefore characterise both $`\chi _k^z`$ and $`\chi _k`$ from a 3D simulation and assess the consequences of using either $`\chi _k^z`$ or $`\chi _k`$ in the calculation of $`P`$.

### 2.3 A Gaussian function for $`\chi _k`$

Stein (1967) and Musielak et al. (1994) suggested several analytical forms for $`\chi _k(\omega )`$. The Gaussian Function (GF hereafter) is the simplest choice and is defined as

$$\chi _k(\omega )=\frac{1}{\omega _k\sqrt{\pi }}e^{(\omega /\omega _k)^2}$$
(8)
where $`\omega _k`$ is its linewidth.

In the time domain, the gaussian function, Eq. (8), is the Fourier transform of a gaussian function whose linewidth is equal to $`2\tau _k`$, where $`\tau _k`$ is a characteristic time correlation length. Hence $`\omega _k`$ and $`\tau _k`$ are related to each other as

$$\omega _k=\frac{2}{\tau _k}.$$
(9)
The characteristic time $`\tau _k`$ is usually associated with the characteristic correlation time-scale of an eddy with wavenumber $`k`$. As in Balmforth (1992), we define it as

$$\tau _k\frac{\lambda }{ku_k}$$
(10)
where the velocity $`u_k`$ of the eddy with wave number $`k`$ is related to the kinetic energy spectrum $`E(k)`$ by (Stein, 1967)

$$u_k^2=_k^{2k}𝑑kE(k).$$
(11)
The parameter $`\lambda `$ in Eq. (10) accounts for our lack of precise knowledge of the time correlation $`\tau _k`$ under stellar conditions.

In the calculation of $`P`$, a GF is usually assumed for $`\chi _k`$ (e.g. Goldreich & Keeley (1977), Balmforth (1992)). This assumption is equivalent to supposing that two distant points in the turbulent medium are uncorrelated.

In Sect. 3, we use a 3D simulation of the upper part of the solar convective zone to derive the $`\nu `$-dependencies of $`\chi _k^z`$ and $`\chi _k`$. Inferred $`\nu `$-dependencies of $`\chi _k^z`$ and $`\chi _k`$ are compared to that of the GF. We next determine several analytical forms for $`\chi _k^z`$ and $`\chi _k`$ that can better represent their $`\nu `$-dependencies.

## 3 Constraints from the 3D simulation

The analysis of a 3D simulation of the upper part of the solar convective zone provides constraints for several physical parameters that enter the theoretical expression for the energy supply rate $`P`$ injected into the solar p modes (Eq. 1). The constraints may be considered to be of two types :

* *static* constraints (*static* refers to spatial and time averages) determine the actual wavenumber dependency of $`E(k,z)`$, the kinetic turbulent spectrum, and $`E_s(k,z)`$, the turbulent spectrum associated with the entropy. The *static* constraints also determine the depth profile of the wavenumber $`k_0^E`$ at which convective energy is injected into the turbulent inertial range of $`E`$ (as in Samadi et al. (2003) we assume that the wavenumber $`k_0^{E_s}`$, at which convective energy is injected into the turbulent inertial range of $`E_s`$, is equal to $`k_0^E`$). They also provide the depth dependence of $`u^2`$, the mean square velocity, $`w^2`$, the mean square vertical component of the velocity, $`\stackrel{~}{s}^2`$, the mean square values of entropy fluctuations and $`\mathrm{\Phi }=u^2/w^2`$, the mean values of the anisotropy studied in Samadi et al. (2003).
* The *dynamic* constraints, on the other hand, concern the frequency component $`\chi _k`$ and $`\chi _k^z`$ (see Sect. 2.2).

The static constraints have been established in Samadi et al. (2003). Here we investigate the dynamic constraints.

### 3.1 The 3D simulation

We study a 3D simulation of the upper part of the solar convection zone obtained with the 3D numerical code developed by Stein & Nordlund (1998).

The simulated domain is 3.2 Mm deep and its surface is 6 x 6 $`\mathrm{Mm}^2`$. The grid of mesh points is 256 x 256 x 163 (i.e. $`23\mathrm{k}\mathrm{m}\times 23\mathrm{k}\mathrm{m}\times 37\mathrm{k}\mathrm{m}`$), the total duration 27 mn and the sampling time 30s.

Outputs of the simulation considered in Samadi et al. (2003) are the velocity field $`u(x,y,z,t)`$ and the entropy $`s(x,y,z,t)`$ where - as in Samadi et al. (2003) \- $`z=rR_{}`$ and $`R_{}`$ is the radius at the photosphere (i.e. where $`T=T_{\mathrm{eff}}`$). The quantities $`u(x,y,z,t)`$ and $`s(x,y,z,t)`$ were used to determine the quantities $`E(k,z)`$, $`E_s(k,z)`$, $`w`$ and $`\stackrel{~}{s}^2`$ involved in the theoretical expression for the excitation rate $`P`$. In the present work we use the velocity field $`u(x,y,z,t)`$ to characterise $`\chi _k`$ and $`\chi _k^z`$.

### 3.2 Fourier transform

At four different layers of the simulated domain, we compute the 3D Fourier transform, with respect to time and in the horizontal plane, of the velocity field $`𝐮`$. These layers cover a region where modes with frequency $`\nu 2`$ mHz are predominantly excited.

This provides $`\widehat{u}(k,z,\nu )`$ where $`k`$ is the wavenumber in the horizontal plane. Next we integrate $`\widehat{u}^2(k,z,\nu )`$ over circles with radius $`k=k`$ at each given layer $`z`$. This yields $`\widehat{u}(k,z,\nu )`$ and therefore $`E(k,\nu ,z)\widehat{u}^2(k,\nu ,z)`$. The quantity $`\chi _k(\nu ,z)`$ is the frequency component of $`E(k,\nu ,z)`$, at fixed $`k`$. Hence, according to Eq.(6) and (7), $`\chi _k(\nu ,z)`$ is obtained from $`E(k,\nu ,z)`$ as:

$$\chi _k(\nu ,z)=\frac{E(k,\nu ,z)}{{\displaystyle d\nu E(k,\nu ,z)}}$$
(12)
where the integration over $`\nu `$ is performed over the frequency range $`[\nu _{\mathrm{max}},\nu _{\mathrm{max}}]`$ corresponding to the window of the Fourier analysis with respect to time ($`\nu _{\mathrm{max}}16`$ mHz).

We proceed in the same manner for the vertical component of $`u`$. This then provides $`E_z(k,\nu ,z)\widehat{u}_{z}^{}{}_{}{}^{2}(k,\nu ,z)`$ and $`\chi _k^z(\nu ,z)`$.

### 3.3 Inferred properties of $`\chi _k`$ and $`\chi _k^z`$

Fig. 1 presents $`\chi _k(\nu )`$ as it is obtained from the simulation for the wavenumber $`k`$ at which $`E(k,z)`$ peaks ($`k=k_0`$). We compare $`\chi _k(\nu )`$ obtained from the 3D simulation analysis with the analytical functions GF, GEF and GLF (Eq. (8), (14), and (15) resp.).

At the top of the superadiabatic region (for instance $`z=0.4`$Mm in Fig.1, this is the layer where the excitation is the largest), the GF does not correctly model $`\chi _k(\nu )`$ (see Fig. 1). However the discrepancies between the GF and the simulation data occur mostly above the solar cut-off frequency ($`\nu 5.5`$ mHz). Discrepancies between the GF and the 3D simulation data have then minor consequences for the p modes excitation in this region. This is not the case deeper in the simulation where the largest discrepancies between the GF and the simulation data occur in the frequency range where the dominant amount of acoustic energy is injected into the p modes ($`\nu 24`$ mHz).

To reproduce the shape of $`\chi _k(\nu )`$ obtained with the 3D simulation, one needs a function which at high frequency decreases more slowly than the GF. For modeling the $`\nu `$-dependency of $`\chi _k(\nu )`$, we thus propose three analytical functions: the Lorentzian function (LF hereafter)

$$\chi _k(\omega )=\frac{1}{\pi \omega _k/2}\frac{1}{1+\left(2\omega /\omega _k\right)^2},$$
(13)
the Gaussian plus an Exponential function (GEF hereafter)

$$\chi _k(\omega )=\frac{1}{2}\left(\frac{1}{\omega _k\sqrt{\pi }}e^{(\omega /\omega _k)^2}+\frac{1}{2\omega _k}e^{|\omega /\omega _k|}\right),$$
(14)
and the Gaussian plus a Lorentzian function (GLF hereafter)

$$\chi _k(\omega )=\frac{1}{2}\left(\frac{1}{\omega _k\sqrt{\pi }}e^{(\omega /\omega _k)^2}+\frac{1}{\pi \omega _k}\frac{1}{1+\left(\omega /\omega _k\right)^2}\right).$$
(15)
All these functions satisfy the condition of normalisation of Eq. (7).

We first assume a constant $`\lambda =1`$. As shown in Fig. 1, all these non-gaussian functions reproduce the $`\nu `$-variation of $`\chi _k`$ better than that obtained using a GF.

In the middle of the excitation region ($`0.5\mathrm{Mm}z0.0\mathrm{Mm}`$) the overall best agreement is obtained with the LF. Below $`z0.5\mathrm{Mm}`$, the LF does not reproduce $`\chi _k`$ well enough but still reproduces its $`\nu `$-variation better than the other models.

However we have so far assumed that $`\lambda `$ (or equivalently the eddy time correlation) is depth independent, which is a strong assumption. When $`\lambda `$ is allowed to vary with $`z`$, we find that decreasing the value of $`\lambda `$ below $`z0.5\mathrm{Mm}`$, the LF best models $`\chi _k`$ below $`z0.5\mathrm{Mm}`$ (e.g. $`\lambda =1.6`$ at z=-0.64 Mm and $`\lambda 1.30`$ at z=-0.99 Mm, see Fig. 2). This shows that the variation with depth of the characteristic time $`\tau _k`$ (or equivalently the characteristic frequency $`\omega _k`$) is not correctly represented by the relations of Eqs. (10-11) when computed assuming a constant $`\lambda `$ below $`z0.5\mathrm{Mm}`$; $`\tau _k`$ increases faster with depth than expected from the relations Eqs. (10-11). It is however found in Sect. 4 that this feature has negligible effect on $`P`$.

The function $`\chi _k^z`$ also decreases with the frequency more slowly than the GF (see Fig. 1). Moreover, decreasing values of $`\lambda `$ for $`z0.5`$ Mm provide a better fit of $`\chi _k^z`$. But in contrast with $`\chi _k`$, $`\chi _k^z`$ is overall better modeled with the GEF for $`z0.5\text{Mm}`$ and with the GLF for $`z0.5\text{Mm}`$ rather than with the LF (not shown).

We conclude from the frequency analysis of the 3D simulation that the simple gaussian function cannot correctly represent the actual dynamic properties of the turbulent medium. One may expect that the GF causes an underestimation of the acoustic energy injected into the solar p modes. Instead the frequency analysis favours a non-gaussian function for $`\chi _k^z`$ and $`\chi _k`$ that decreases more slowly with $`\nu `$ than the GF.

## 4 Consequences in terms of p modes excitation

### 4.1 Computations of the excitation rate $`P`$

Computation of the excitation rate $`P`$ is performed as in Samadi et al. (2003) except that here two analytical functions other than the GF are assumed for $`\chi _k`$, as discussed in Sect. 3. The computation process is summarised as follows: The eigenfunctions ($`\xi _\mathrm{r}`$) and their frequencies ($`\nu `$) are computed with Balmforth’s (1992) non-adiabatic code for a solar 1D mixing-length model based on Gough’s (1977) non-local time-dependent formulation of convection.

The quantities $`\mathrm{\Phi }`$, $`w^2`$ and $`s^2`$ are obtained from the 3D simulation. The $`k`$-dependency of $`E(k,z)`$ is the Extended Kolmogorov Spectrum (EKS hereafter) defined as :

$$\begin{array}{cccc}E(k)(k/k_0)^{+1}\hfill & \text{for}& k_0>k>k_{\mathrm{min}}& \\ E(k)(k/k_0)^{5/3}\hfill & \text{for}& k>k_0& \end{array}$$
(16)
In Eq. (16), the wavenumber $`k_0`$ is the wavenumber at which $`E(k)`$ peaks and $`k_{\mathrm{min}}`$ is the minimal wavenumber reached by the 3D simulation ($`k_{\mathrm{min}}=1.05\mathrm{Mm}^1`$). The variation with depth of $`k_0`$ is also given by the 3D simulation. The $`k`$-dependency of the EKS reproduces the global features of $`E`$ arising from the 3D simulation. The same model is considered for $`E_s(k,z)`$. $`E(k,z)`$ and $`E_s(k,z)`$ satisfy the normalisation conditions:

$$\begin{array}{ccc}_{k_{\mathrm{min}}}^+\mathrm{}dkE(k,z)\hfill & =\hfill & 1/2\mathrm{\Phi }w^2\hfill \\ _{k_{\mathrm{min}}}^+\mathrm{}dkE_s(k,z)\hfill & =\hfill & 1/2\stackrel{~}{s}^2.\hfill \end{array}$$
(17)
The total energy contained in $`E(k,z)`$ and $`E_s(k,z)`$, and their depth dependencies, are then obtained from the 3D simulation according to Eq.(17). These theoretical estimates for $`P`$ are then compared with the ‘observed’ $`P`$ from Chaplin et al. (1998)’s seismic data, calculated according to the relation :

$$P(\omega _0)=2\eta \frac{I}{\xi _\mathrm{r}^2(r_s)}v_s^2(\omega _0)$$
(18)
where $`r_s`$ is the radius at which oscillations are measured,

$$I_0^M𝑑m\xi _\mathrm{r}^2$$
(19)
is the mode inertia and where the mode damping rate ($`\eta `$) and the mode surface velocity ($`v_s`$) are obtained from Chaplin et al. (1998). In Eq. (18), the mode mass $`I/\xi _\mathrm{r}^2(r_s)`$ is given by the GMLT model and we adopt $`r_s=R_{}+200`$ km consistently with Chaplin et al. (1998)’s observations.

### 4.2 Comparisons with observations

We investigate the effect of using different analytical functions for $`\chi _k`$ (Sect. 3) in the computation of $`P`$ . We first assume a constant value $`\lambda =1`$. Results are shown in Fig. 3. Computations performed with the GF underestimate the observed $`P`$ values by a factor $`2.7`$. On the other hand, the LF, GEF and GLF choices result in larger values for the computed $`P`$ than the GF one ($`2`$ times larger). This brings them closer to the observations, compared with the GF choice for $`\chi _k`$. The reason is that all the non-gaussian functions (the LF , the GEF and the GLF) \- which indeed better model $`\chi _k^z`$ and $`\chi _k`$ from the 3D simulation than does the GF - decrease more slowly with $`\nu `$ than the GF in the frequency range where the mode amplitudes are large ($`\nu 24`$ mHz). Consequently a larger amount of acoustic energy is injected into the modes with the non-gaussian functions than with the GF.

In Sect. 3.3, the overall best models for $`\chi _k`$ were obtained with the LF and with decreasing values of $`\lambda `$ below $`z0.5`$ Mm. We use a simple model for the depth variation of $`\lambda `$:

$$\begin{array}{cccc}\lambda \hfill & =1\hfill & \text{for}\hfill & z>0.5\mathrm{Mm}\hfill \\ \lambda \hfill & =0.9+0.71(0.64\text{Mm}+z)\hfill & \text{for}\hfill & 0.5z1\mathrm{Mm}\hfill \\ \lambda \hfill & =0.35\hfill & \text{for}\hfill & z<1\mathrm{Mm}\hfill \end{array}$$
(20)
We have computed $`P`$ according to the simple model of Eq. (20) and assuming the LF. We find no significant changes for $`P`$ compared to the calculations in which a constant value $`\lambda =1`$ is assumed (not shown).

In Sect. 3 we found that $`\chi _k^z`$ is better modelled with the GEF for $`z0.5\text{Mm}`$ and with the GLF for $`z0.5\text{Mm}`$ rather than with the LF. However, as the stochastic excitation is the largest in the range $`0.5\text{Mm}z0`$ Mm, we can assume the GEF in all the domain. The LF results in a value for $`P_{\mathrm{max}}`$ slightly larger than the one resulting from the GEF (only $`1.2`$ larger). A better agreement is then obtained with the analytical functions which fits best $`\chi _k`$ (i.e. the LF) than the one which fits best $`\chi _k^z`$ (i.e. the GEF) in contrast with the intuitive idea mentioned in Sect. 2.2 that the excitation of radial p modes depends rather on the properties of $`\chi _k^z`$ than on those of $`\chi _k`$.

In the frequency range where observational constraints are available, differences between results obtained with the different adopted non-gaussian functions are of the same order as the actual error bars associated with the observed $`P`$ values.

Below the frequency range of the observations - i.e. below $`\nu 1.8`$mHz - the differences between the different non-gaussian functions are very large compared to the error bars (see bottom panel of Fig. 3). Those differences are related to differences in the $`\nu `$-variation of the non-gaussian models. Observational constraints at low frequency could therefore confirm that the LF is indeed the best representation for $`\chi _k(\nu )`$.

Important discrepancies still remain at high frequency ($`\nu 3.5`$ mHz). The excitation rate derived from the observations decreases as $`\nu ^{6.2}`$ above $`\nu 3.5`$ mHz) whereas the computed $`P`$ decreases as $`\nu ^1`$ (see bottom panel of Fig. 3). Possible origins of this discrepancy are discussed in Sect. 5.5.

Another consequence of a non-gaussian dependence of $`\chi _k`$ (or $`\chi _k^z`$) with the frequency is a larger relative contribution of the Reynolds stress $`P_R`$ to the mode excitation rate $`P`$. This is shown in Fig. 4. The GF generates a relative contribution of the Reynolds stress to the excitation ($`P_R/P`$) which is smaller than that obtained assuming a non-gaussian function (e.g. for the LF, the relative contribution of $`P_R`$ to the excitation is at least $`2`$ times larger than with the GF). Excitation by the entropy fluctuations takes place predominantly at the top of the excitation region over a thin layer ($`0.2`$ Mm) while that due to the Reynolds stress extends deeper below ($`0.52`$ Mm). At the top of the excitation region, the discrepancy between the GF and $`\chi _k`$ inferred from the 3D simulation mainly occurs above $`\nu 5`$ mHz and thus has a small impact on mode excitation. This is not the case deeper in the excitation region where the GF under-estimates $`\chi _k`$ in a frequency range increasing inward.

## 5 Conclusion and discussion

### 5.1 A Non-Gaussian eddy time correlation

In the present work we characterize *empirically* $`\chi _k`$ and of $`\chi _k^z`$, the frequency components of the correlation product of the turbulent velocity field and of its vertical component respectively. A frequency analysis of a solar 3D simulation shows that at large scales ($`kk_0`$) the gaussian function significantly underestimates $`\chi _k`$ and $`\chi _k^z`$ in the frequency range ($`\nu 24`$ mHz) where acoustic energy injected into the solar p modes is the largest.

As a result, the maximum value of $`P`$ is found $`2.7`$ smaller than the solar seismic constraints.

This partly explains the *underestimate* of the values of solar p mode excitation rates obtained by Houdek et al. (1999) whose computations are based on the theoretical expression by Balmforth (1992) and by Samadi et al. (2003).

In order to reproduce the main properties of $`\chi _k`$ (or $`\chi _k^z`$) , one has to consider a model which must decrease more slowly with $`\nu `$ than the GF. We then assume for $`\chi _k`$ and $`\chi _k^z`$ three different simple analytical forms : the Lorentzian Function (LF), the so-called ’Gaussian Exponential Function’ (GEF, which is composed by the GF plus an exponential function) and the so-called ’Gaussian Lorentzian Function’ (GLF, which is composed by the GF plus a Lorentzian function).

From the top of the excitation region (which corresponds to the top of the superadiabatic region) down to the middle of the excitation region ($`z0.5`$ Mm where $`z`$ is the distance to the radius at the photosphere), the best agreement between $`\chi _k`$ and the analytical approximations is obtained with the LF and with $`\lambda =1`$. Deeper within the excitation region ($`z0.5`$Mm), the agreement is better for $`\lambda <1`$ and $`\lambda `$ decreasing with depth.

The frequency dependencies of $`\chi _k^z`$ and $`\chi _k`$ are found to be very similar. However $`\chi _k^z`$ is best modeled by the GEF. As for $`\chi _k`$, the agreement is better below $`z0.5`$ Mm with decreasing values of the parameter $`\lambda `$ than with $`\lambda =1`$.

Assuming a non-gaussian function - either the LF, the GEF or the GLF - results in values for $`P_{\mathrm{max}}`$, the maximum of excitation power, which are $`2`$ times larger than when assuming the GF and brings $`P_{\mathrm{max}}`$ much closer to the maximum of $`P`$ derived from the solar seismic data of Chaplin et al. (1998).

We also find that taking into account the variation of $`\lambda `$ with depth for $`z`$ below $`0.5`$ Mm does not significantly change the values of $`P`$. A constant value can then be assumed in the calculation of the solar p mode excitation rates. The constant value of $`\lambda `$ on the other hand plays an important role and we find $`\lambda =1`$.

We have investigated the sensitivity to the adopted representation for $`\chi _k`$ : Although the LF fits best the $`\nu `$-variation of $`\chi _k`$ inferred from the 3D simulation, the GLF results in value for $`P_{\mathrm{max}}`$ closer to the seismic constraints. However, the differences obtained with the different non-gaussian approximations for $`\chi _k`$ are globally smaller than the actual error bars associated with the observations of Chaplin et al. (1998). On the other hand, below the frequency range where observational constraints on $`P`$ are available (i.e. below $`\nu 1.8`$mHz), the differences between $`P`$ obtained with different non-gaussian functions are very large compared to the current error bars. Those differences are directly related to the diffences in the $`\nu `$-variation of the non-gaussian forms investigated in this work. This suggests that accurate enough data below this frequency range, could provide confirmation that the LF is indeed the best model for $`\chi _k`$.

### 5.2 Relative contribution of the entropy fluctuations to the excitation

The non-gaussian character of $`\chi _k`$ causes the excitation region to extend deeper ($`500`$ km for modes of order $`n=20`$) than with the GF ($`200`$ km resp.). The largest entropy fluctuations mainly occur at the outermost part of the convective zone (CZ) over a very thin region ($`100`$ km) while excitation by the Reynolds stress contribution occurs on a more extended region. Consequently the non-gaussian property of $`\chi _k`$ leads to a relatively larger contribution of the Reynolds stress to the excitation than in the case of a GF. As a result, the Reynolds stress contribution is of the same order as the contribution arising from the advection of the turbulent fluctuations of entropy by the turbulent movements (the so-called entropy source term). This is in contrast with previous results (Samadi et al., 2001) based on the GF which concluded that the entropy source term dominates the Reynolds stress by about $`20`$. It also differs with results by Goldreich et al. (1994) who found that the excitation arising directly from the entropy fluctuations dominates by about $`10`$.

On the other hand, in Stein & Nordlund (2001), the excitation by turbulent pressure (Reynolds stress) is found dominant ($`4`$ times larger) whereas here we find that the contribution of the entropy source term cannot be neglected. Whether this is the signature of some deficiency in the present excitation model is an open question.

### 5.3 Summary

We show that the usually adopted *gaussian function* for $`\chi _k`$ is neither consistent with the properties of $`\chi _k`$ inferred from the 3D simulation nor does it reproduce the observed maximum of the solar p-modes excitation rates.

Following an empirical approach we improve the model of the convective eddy time-correlation $`\chi _k`$ which enters the current model of stochastic excitation. We then show that to reproduce *both* the $`\nu `$-variation of $`\chi _k`$ as inferred from the 3D simulation and the observed maximum of the solar p-modes excitation rates one has to consider a non-gaussian form which decreases at high frequency slower than the GF, as do the different non-gaussian functions investigated here.

The use of non-gaussian functions, for instance the LF, reproduces reasonably well the maximum value of the rate at which solar p-modes are excited *without any ajdustements of free parameters* or *without introducing a scaling factor*, in contrast with previous approaches (e.g. Balmforth, 1992; Goldreich et al., 1994; Samadi et al., 2001). We then solve the problem of the underestimation by the previous theoretical approaches. Furthermore the use of such a non-gaussian form for $`\chi _k`$ makes the contribution of the turbulent pressure to the excitation much larger than in previous works making our results more consistent with that by Stein & Nordlund (2001).

Our investigation clearly emphasizes the *non-gaussian character of the solar p modes excitation* as a result of the *non-gaussian property of the convective eddies time-correlations*. It also shows that the dynamic properties of the solar turbulent convection inferred from the 3D simulation are consistent with the helioseismic data.

We stress that only simple non-gaussian forms for $`\chi _k`$ have been investigated here. More sophisticated forms are likely to improve the agreement with the $`\nu `$-dependency of $`\chi _k`$ (or $`\chi _k^z`$). This would not affect the main conclusions presented in the present paper.

### 5.4 Possible origin of the non-gaussian property of $`\chi _k`$

We recall that $`\chi _k`$ measures the temporal evolution of the correlation between two points of the turbulent medium separated by a distance of $`2\pi /k`$. A gaussian time-correlation means that the fluid motions in the medium are random in time. Departure from a gaussian time-correlation at large scales ($`kk_0`$) suggests that a strong correlation exists at that scale.

Downward plumes are likely to be responsible for the non-gaussian behaviour of $`\chi _k`$. Downward and upward convective motions are indeed highly asymmetric (Stein & Nordlund, 1998): downward flows are associated with patterns (plumes) which are more coherent than the upward moving structures (Rieutord & Zahn, 1995). The upward flows are associated with less coherent and more random structures (granules) characterised by a broad variety of sizes and lifetimes (Rieutord & Zahn, 1995). The non-gaussian behavior of $`\chi _k`$ can most probably be attributed to plumes. This however remains to be checked (work in progress).

### 5.5 Possible origin of the remaining discrepancy

Despite a clear improvement in the agreement between observed and theoretical excitation rates, important discrepancies between the computed $`P`$ and the solar measurements still remain at high frequency $`\nu 3.5`$ mHz (see Sect. 4.2 and Fig. 3).

On the ‘observational side’, at high frequency, larger uncertainties for the damping rates $`\eta `$ induce larger uncertainties on the derived supply energy rates.

On the theoretical side, part of the discrepancy might well be attributed to a poor description of the eigenfunctions at high frequency. Indeed, the discrepancies between the calculated eigenfrequencies and the observed ones are largest at high frequency ($`\nu 3`$mHz). This indicates that the description of the eigenfunctions are less accurate at high frequency. As the expression for calculating $`P`$ involves the first and second derivatives of the mode eigenfunction, the lack of accuracy in the calculation of the eigenfunctions has a larger impact on $`P`$ at high frequencies than at small frequencies.

Other possible causes can perhaps be related to our simplified excitation model which assumes isotropic turbulence. Indeed the current theory assumes that the stochastic excitation is the same in all three directions, particularly between the ascending and descending flows. However the kinetic energy and entropy fluctuations are larger in the downward flows than in the upward flows (Stein & Nordlund, 1998). Therefore the driving arising from the advection of the turbulent fluctuations of entropy by the turbulent movements differs significantly between the elements moving downwards and those moving upwards. As the entropy fluctuations are largest in the outermost part of the convective zone, the above mentioned asymmetry will predominantly affect the high frequency modes.

Moreover, it is also assumed that the total kinetic energy, $`E`$, is isotropically injected in all 3 directions. Excitation of the radial p modes results from the vertical component of the velocity. However at the top of the convective zone, the distribution of kinetic energy in $`E(k,z)`$ and in $`E_z(k,z)`$ are very different from each other. These differences may affect more strongly the high frequency modes. Consequences of these departures from the isotropic assumption need to be further investigated.

### 5.6 Perspectives

The non-gaussian property of $`\chi _k`$ and its consequences for the stochastic excitation has been investigated so far only for the *Sun*. However such a non-gaussian feature of the turbulence will most likely also be of importance for solar-like oscillating stars more massive than the Sun, provided our analysis is also valid for these stars. This can substantially change the excitation spectrum $`P`$ for such stars compared to that which is currently predicted.

Therefore investigations of p mode excitation in hotter and more massive stars must be undertaken, which should proceed in two steps : first, the validity of the present results obtained in the solar case must be investigated for other stars with, for instance, the help of dedicated 3D simulations. The conclusions which will drawn from this first step must be used in a second step to study the frequency dependence and the magnitude of $`P`$ for different solar-like oscillating stars (see preliminary results in Samadi et al., 2002).

Future space missions such as COROT (Baglin & The Corot Team, 1998), MOST (Matthews, 1998) and Eddington (Favata et al., 2000) will provide high-quality data on seismic observations. COROT will be the first mission that will provide high precision mode amplitudes and linewidths in other stars. This high-quality data will allow us to derive the excitation rate $`P`$. and will provide improved observational constraints on the theory of stochastic excitation which is, at present, poorly constrained by observation.

###### Acknowledgements.

We thank H.-G. Ludwig for useful discussions and valuable help in analysing the simulated data and we are indebted to G. Houdek for providing us the 1D solar model. We gratefully thank our referee (T. Appourchaux) for pertinent remarks which led to an improvement of the paper. RS’s work has been supported in part by the Particle Physics and Astronomy Research Council of the UK under grant PPA/G/O/1998/00576.