Numerical simulations of the kappa-mechanism with convection

0.46 MB·English

by  T. Gastine

AstrophysicsandSpaceScience
DOI10.1007/s•••••-•••-••••-•
Numerical simulations of the κ-mechanism with
convection
T. Gastine1 • B. Dintrans1
0
1
0
2

n
a
J (cid:13)c Springer-Verlag••••

5
Abstract A strong coupling between convection and These driving zones should nevertheless be located at

] pulsations is known to play a major role in the dis- a precise radius inside the star, neither too close to
R
appearance of unstable modes close to the red edge of the surface nor to deep, in order to balance the damp-
S
the classical Cepheid instability strip. As mean-field ing that mainly occurs at the surface. It defines the
.
h models of time-dependent convection rely on weakly- so-called transition region that shapes the limit be-
p
constrainedparameters(e.g.Baker1987),wetacklethis tweenthequasi-adiabaticinteriorandthestronglynon-
-
o problem by the means of 2-D Direct Numerical Simu- adiabatic surface. For classical Cepheids that pul-
tr lations (DNS) of κ-mechanism with convection. sate in the fundamental acoustic mode, this transi-
s Using a linear stability analysis, we first determine tion region is located at a temperature T ≃ 4 ×
a
[ the physical conditions favourable to the κ-mechanism 104 K corresponding to the second helium ionisation
to occur inside a purely-radiative layer. Both the in- (Baker and Kippenhahn 1965). However, the bump lo-
2
stabilitystripsandthenonlinearsaturationofunstable
v cationisnotsolelyresponsiblefortheacousticinstabil-
modes are then confirmed by the corresponding DNS.
7 ity. A careful treatment of the κ-mechanism would in-
8 We next present the new simulations with convection, volve dynamical couplings with convection, metallicity
4 where a convective zone and the driving region over-
effects and both realistic equations of state and opac-
0
lap. The coupling between the convective motions and
. itytables(Bono, Marconi, and Stellingwerf1999). The
1 acoustic modes is then addressed by using projections
purpose of our model is to simplify the hydrodynamic
0 onto an acoustic subspace.
0 approachwhileretainingtheleading-order phenomenon
1 –e.g. thelocationoftheopacitybumporitsamplitude–
Keywords Hydrodynamics - Instabilities - Stars: os-
v: cillations - Convection - Methods: numerical suchthatfeasibledirectnumericalsimulationsoftheκ-
i mechanism with convection can be achieved.
X
r
a 1 Introduction
2 The first step: radiative models in 1-D
Eddington (1917) discovered an excitation mechanism
of stellar oscillations that is related to the opacity be-
Wefirstfocusonradialmodespropagatinginapurely-
haviour in ionisation zones: the κ-mechanism, where
radiative and partially-ionised layer and then restrict
κ denotes the opacity. This mechanism occurs when
our study to the 1-D case. Our model consists of a
the opacity varies during compression phases so as
local zoom about an ionisation region and is composed
to block the emerging radiative flux (Zhevakin 1963;
byamonatomicandperfectgas(withγ =c /c =5/3)
p v
Cox and Whitney1958). Ionisationregionscorrespond
underbothaconstantgravity~gandkinematicviscosity
to a strong increase in opacity that leads to the opac-
ν. All quantities in our model are dimensionless, e.g.
ity bumps responsibleforthelocalexcitationofmodes.
temperatureisgiveninunitofthesurfacetemperature
ofthestar,lengthinunitofafractionofthestarradius
T.Gastine and velocity in unit of the sound speed at the surface.
B.Dintrans In order to easily investigate the influence of the
1LATT,Universit´edeToulouse, CNRS,14 avenueEdouard Be- opacity bump on the stability, we adopt the following
lin,F-31400Toulouse,France
2
Fig. 1 Influence of the hollow parameters on the conduc-
−2
tivity profile for Kmax =10 and Tbump/Tsurf =3.5.
conductivity hollow to mimic this bump (as κ ∝ 1/K,
see Gastine and Dintrans 2008a, hereafter Paper I): Fig. 2 a) Vertical profile of the Ψ coefficient (Eq. 3).
Thesuperimposedgreenzonerepresentsthelocationwhere
Ψ=1±0.5. b) The real part of the work integral plotted
−π/2+arctan(σT+T−) for the three different equilibrium models. The dot-dashed
K (T )=K 1+A , (1)
0 0 max(cid:20) π/2+arctan(σe2) (cid:21) black line is for the constant radiative conductivity case.
c) Corresponding radiative conductivity profiles. d) Corre-
withA=(Kmax−Kmin)/Kmax andT± =T0−Tbump± sponding equilibrium fields dKT/dz (Eq. 2).
e. HereT istheequilibriumtemperatureprofile,T
0 bump
is the hollow central temperature, and σ, e and A de-
luminosity. Ψ is the ratio between the thermal energy
note its slope, width and relative amplitude, respec-
embeddedbetweenthegivenradiusandthesurfaceand
tively. Examplesofcommonvaluesoftheseparameters
the energy radiated during an oscillation period.
are provided in Fig. 1.
2.2 The linear stability analysis
2.1 Conditions for the instability
Inordertocheckthesetwoconditions,weperformalin-
Using the well-known work integral formalism (e.g.
earstabilityanalysisbyusingtheLinearSolverBuilder
Unno et al. 1989), it is easy to show that the follow-
code(Valdettaro et al.2007). Bystartingwiththecon-
ing condition is necessary to get unstable modes by κ-
ductivity profiles shown in Fig. 1, we first compute the
mechanism in variable stars:
equilibrium setup that is solution of both the hydro-
static and radiative equilibria given by
dK ∂lnK
T 0
<0 where K ≡ . (2)
T
dz (cid:18)∂lnT (cid:19)
0 ρ0 ∇~p =ρ ~g and div[K (T )∇~T ]=0, (4)
0 0 0 0 0
However, this condition is not sufficient as the ionisa-
tion region should also be located at a given location where p0 and ρ0 denote the equilibrium pressure and
in the star, neither too deep nor too close to its sur- density. Then, we solve the linear oscillation equations
face. This area is the so-called transition region that for the perturbations by seeking normal modes of the
separatesthequasi-adiabaticinteriorfromthestrongly form exp(λt) with λ = τ +iω, where ω is the mode
non-adiabatic surface and is defined by (e.g. Cox 1980) frequency and τ its growth or damping rate (unstable
modes corresponding to τ > 0). By varying the hol-
low parameters, we then determine the physical con-
hc T i∆m
Ψ= v 0 ≃1, (3) ditions that lead to unstable modes excited by the κ-
ΠL
mechanism in our layer (cf. Paper I).
where ∆m is the integrated mass between the consid- The two criteria (2) and (3) predict that modes are
eredpointandthesurface,ΠthemodeperiodandLthe unstable for a “sufficient” hollow (adequate amplitude
Numericalsimulationsoftheκ-mechanismwithconvection 3
and width) that is located in the transition region.
Using the work integral computed from the obtained
eigenvectors,weindeedrecoverinstabilitystripsforthe
κ-mechanismbyconsideringthethreedifferentequilib-
rium setups (Fig. 2):
• The first case (dotted blue line, T /T = 1.7)
bump surf
corresponds to a hot star with an ionisation region
close to the surface. As the surrounding density is
small, Ψ≪1 and the conductivity hollow hardly in-
fluencestheworkintegral: drivingisunabletoprevail
over damping, leading to τ <0 or stable modes.
• In the second case (solid green line, T /T =
bump surf
2.1), the radiative conductivity begins to decrease
significantly at the location of the transition region
whereΨ≃1. Asaconsequence,drivingbecomesim-
portant. Furthermore, the radiative conductivity in-
creaseoccurswherenon-adiabatic effectsarealready
significant, i.e. Ψ<1there. Itmeansthatnodamp-
ing occurs between the hollow position and the sur-
face astheradiative fluxperturbations arefrozen-in.
Drivingisovercomingdampinginthiscase,therefore
τ >0 and modes are unstable.
Fig. 3 a) Temporal power spectrum for the momentum
• The third case (dashed red line, T /T = 2.8)
bump surf in the (z, ω)-plane. b) The resulting spectrum after in-
corresponds to a cold star of which ionisation re- tegrating over depth. c) Comparison between normalised
gionislocateddeeperinthestellaratmospherewhere momentum profiles for n = (0, 2) modes according to the
Ψ ≫ 1. Ionisation then occurs in a quasi-adiabatic DNS power spectrum (solid black line) and to the linear
location. As a consequence, the excitation provided stability analysis (dotted blue line).
bytheconductivityhollowcannotbalancethedamp-
ing arising at the top of the layer, thus τ < 0 and
To determine which modes are present in the DNS
modes are stable.
in the nonlinear-saturation regime, we first perform
In conclusion, our first radiative model of the κ- a temporal Fourier transform of the momentum field
mechanism is well able to reproduce the main physics ρu(z,t) and plot the resulting power spectrum in the
involved in the oscillations of Cepheid stars: (i) one (z,ω)-plane (Fig. 3a). The mean spectrum follows by
should have dK /dz < 0 to drive the oscillations; (ii) integrating ρu(z,ω) over depth (Fig. 3b). With this
T
the thermal engine underlying the conductivity hollow method, acoustic modes are extracted because they
c
should also be located in the transition region where emerge as “shark-fin profiles” about definite eigenfre-
Ψ≃1 (Fig. 2). quencies (Dintrans and Brandenburg 2004).
Several discrete peaks corresponding to normal
2.3 Study of the nonlinear saturation modes well appear but the fundamental mode close
to ω = 5.439 clearly dominates. Moreover, the lin-
0
To confirm both the results obtained previously in the
ear eigenfunctions are compared to the mean profiles
linear stability analysis and to investigate the result-
computed from a zoom taken in the DNS power spec-
ing nonlinear saturation, we perform direct numerical
trumabouteigenfrequenciesω =5.439andω =11.06
0 2
simulations of our problem. The idea is to start from
(Fig.3c). Theagreementbetweenthelineareigenfunc-
the most favourable initial conditions found in the sta-
tions (dotted blue lines) and the DNS profiles (solid
bility analysis and then advance the general hydro-
black lines) is remarkable. In summary, Fig. 3 shows
dynamic equations in time. We use a modified ver-
that several overtones are present in the DNS, even for
sionofthepublic-domainfinite-differencePencilCode2
long times. These overtones are however linearly stable
that includes an implicit solver of our own for the
suggestingthatsomeunderlyingenergytransfersoccur
radiative diffusion term in the temperature equation
between modes through nonlinear couplings.
(Gastine and Dintrans 2008b, hereafter Paper II).
These nonlinear interactions are investigated thanks
toatoolthatmeasuresthesoundgenerationbyturbu-
2Seehttp://pencil-code.googlecode.com/. lent flows (e.g. Bogdan, Cattaneo, and Malagoli 1993).
4
It involves projections of the DNS physical fields onto 3 Convective models in 2-D
the regular and adjoint eigenvectors that are solutions
to the linear oscillation equations. The projection co- 3.1 From radiative to convective setups
efficients then give the time evolution of each acoustic
After the previous results obtained on purely radia-
mode that is present in the DNS separately. Further-
tive models in 1-D, we next address the convection-
more, as the evolution of the kinetic energy of each
pulsation coupling that is suspected to quench the ra-
mode is also accessible with this method, one can fol-
dial oscillations of Cepheids close to the red edge. The
low the energy transfer between modes.
idea is to slightly modify the unstable radiative setup
to obtain a convective zone superimposed to the ion-
isation region in 2-D simulations. The convective in-
stability obeys Schwarzschild’s criterion given by (e.g.
Chandrasekhar 1961)
dT g gK (T )
0 0 0
> =⇒ F > . (5)
bot
(cid:12) dz (cid:12) c c
(cid:12) (cid:12) p p
(cid:12) (cid:12)
(cid:12)where(cid:12) F = K |dT /dz| is the imposed bottom flux.
bot 0 0
For a given gravity g and conductivity hollow K (T ),
0 0
convectionwillthereforedevelopforalargeenoughbot-
tom flux satisfying Eq.(5).
Fig. 4 Time evolution of the kinetic energy ratio for the
n=0 and n=2 acoustic modes.
Figure 4 shows an example of the time evolution
of the kinetic energy ratio E /Etot for the two modes
n kin
n = 0 (i.e. the fundamental one) and n = 2 (i.e. the
second overtone). Here E denotes the kinetic energy
n
embeddedintheacousticmodeofordern,whileEtot is
kin
the total kinetic energy in the simulation. After its lin-
ear transient growth, agiven fraction of the fundamen-
tal mode energy is progressively transferred to the sec-
ond overtone until the nonlinear saturation is achieved
abovetimet≃150,withfinally85%ofthetotalkinetic
Fig.5 Snapshotofthevorticityfield∇~ ×~vina2-Dmodel
energyinthen=0modeandtheremaining15%inthe
with convection and κ-mechanism.
n=2 one. We recall that it is at a first glance surpris-
ingtodetectthissecondovertoneatlongtimesbecause Figure5displaysasnapshotofthevorticityfieldob-
it is linearly stable and should normally be damped by tainedinsuchaconvectivesimulationattimet=2500.
diffusive effects. In order to ensure that the thermal relaxation is well
Thereasonforitspresenceliesinafavored coupling achieved,weranthissimulationforamuchlongertime
in the period ratio between these two modes. Indeed, than in the purely-radiative case and the resolution is
the fundamental period is P0 = 2π/ω0 ≃ 1.155 in this 256×256 gridpoints. Asexpected, convection develops
simulation while the n = 2 one is P2 ≃ 0.568. As a in the middle of the box where the conductivity hollow
consequence, the corresponding period ratio is close to is located1. Moreover, because of the density contrast
one half (P2/P0 ≃ 0.491) such that the second over- across the convection zone, the vorticity is trapped in
tone is involved in the nonlinear saturation through a strong convective plume downdrafts that easily over-
2:1 resonance with the fundamental mode. The lin- shoot in the stable zone below due to their low P´eclet
ear growth of the fundamental mode is then balanced number (Dintrans 2009).
by the pumping of energy to the stable second over-
tone that behaves as an energy sink, leading to the full
1An animation of this simulation is provided at
limit-cycle stability.
http://www.ast.obs-mip.fr/users/tgastine/conv_kappa.avi.
Numericalsimulationsoftheκ-mechanismwithconvection 5
3.2 Evolution of acoustic modes with convection (dashed green line). Indeed, after a long transient dur-
ing which it reaches a non-trifling value (∼ 30%), it
To address the time evolution of the acoustic modes progressively tends to zero at long times while the con-
thatarepresentintheDNSwithconvection,weproject vective motions are responsible for the whole kinetic
as before the physical fields onto an appropriate acous- energy content. The radial oscillations excited by κ-
tic subspace to get the projection coefficient c (t). mechanism are thus quenched by convection and this
ℓn
Here ℓ denotes the mode degree while n is its order situation is relevant to the red edge where the unstable
(with k =(2π/L )ℓ and L is the width of the box). acoustic modes are known to be damped by convective
x x x
motions occuring at the star’s surface.
We thus have shown that the Rayleigh number is
a key parameter for this problem. But further inves-
tigations of the convection-pulsation coupling are in
progress, by especially checking time-dependent theo-
riesofconvection forwhichseveralfreeparameterscan
be constrained by such kind of direct numerical simu-
lations (e.g. Yecko, Kolla´th, and Buchler 1998).
Acknowledgements Calculations were carried out
on the CALMIP machine of the University of Toulouse
(http://www.calmip.cict.fr/).
References
Fig. 6 Evolution of the acoustic kinetic energy ratio for Baker,N.,Kippenhahn,R.: Astrophys.J.142,868(1965)
twosimulationswithaweakconvection(solidblueline)and Baker, N.H.: In: Hillebrandt, W., Meyer-Hofmeister, E.,
a stronger one (dashed green line). Thomas,H.C.(eds.)PhysicalProcessesinComets,Stars
and Active Galaxies, p. 105 (1987)
Bogdan, T.J., Cattaneo, F., Malagoli, A.: Astrophys. J.
Fig. 6 shows the total energy embedded in acoustic
407, 316 (1993)
modes for two simulations with convection that differ Bono, G., Marconi, M., Stellingwerf, R.F.: Astrophys. J.
by their Rayleigh number given by Suppl. Ser. 122, 167 (1999)
Chandrasekhar, S.: Hydrodynamic and hydromagnetic
stability. Oxford University Press: Clarendon, Oxford
Ra=−gd4conv ds = (∇−∇adia)gd4conv, (6) (1961)
νχc dz νχH Cox,J.P.: Theoryofstellarpulsation.PrincetonUniversity
p p
Press, Princeton (1980)
where H = −(dlnp /dz)−1 is the pressure scale Cox, J.P., Whitney, C.: Astrophys. J. 127, 561 (1958)
p 0
Dintrans,B.: CommunicationsinAsteroseismology158,45
height, d the vertical extent of the convection zone,
conv (2009)
χ=K /ρ c theradiativediffusivityandstheentropy.
0 0 p Dintrans, B., Brandenburg, A.: A&A 421, 775 (2004)
Thisnumberallowstoquantifythestrengthofthecon- Eddington, A.S.: The Observatory 40, 290 (1917)
vective motions, that is, the critical Rayleigh number Gastine, T., Dintrans, B.: Astron. Astrophys. 484, 29
above which convection develops is about 103 for poly- (2008a)
Gastine, T., Dintrans, B.: Astron. Astrophys. 490, 743
tropic stratifications (e.g. Gough et al. 1976).
(2008b)
Inthefirstsimulation(Ra=4×104,solidblueline), Gough, D.O., Moore, D.R., Spiegel, E.A., Weiss, N.O.: As-
the Rayleigh number is slightly super-critical and the trophys. J. 206, 536 (1976)
convection intensity is weak. In this case, we see that Unno, W., Osaki, Y., Ando, H., Saio, H., Shibahashi, H.:
the acoustic energy ratio remains large (i.e. & 80%, Nonradialoscillationsofstars.UniversityofTokyoPress,
Tokyo (1989)
the remaining 20% being in the convection) and the
Valdettaro, L., Rieutord, M., Braconnier, T., Fraysse, V.:
nonlinear saturation is reached in the same way than
JCoAM 205(1), 382 (2007)
in the previous radiative models. In other words, the Yecko, P.A., Kolla´th, Z., Buchler, J.R.: Astron. Astrophys.
modepropagationisnotaffectedbyconvectivemotions 336, 553 (1998)
and the convection can be seen as frozen-in. Zhevakin, S.A.: Annu. Rev. Astron. Astrophys. 1, 367
(1963)
On the contrary, increasing the Rayleigh number
leads to a more vigorous convection and the acoustic
energyratioismuchmoreaffectedbyconvectiveplumes ThismanuscriptwaspreparedwiththeAASLATEXmacrosv5.2.