A&A 406, 657-666 (2003)
DOI: 10.1051/0004-6361:20030739
C. Weidner - D. Koester
Institut für Theoretische Physik und Astrophysik, Universität Kiel, 24098 Kiel, Germany
Received 7 March 2003 / Accepted 28 April 2003
Abstract
The numerical simulation method for variable white
dwarfs of Ising & Koester (2001) is extended to variable DB stars
(DBVs or V777 Herculis stars). We find the same general behavior of
nonlinear effects, with a sinusoidal light-curve at small amplitudes,
strongly non-sinusoidal variations at intermediate amplitudes and a
change back to smaller and sinusoidal variations at the largest
driving pressure amplitudes. The transition between the various
regimes is however shifted to higher amplitudes compared to DAVs
because of the smoother reaction of the DB convection zones to
perturbations. A peculiar event of the prototype GD358 - with all
pulsation power going into one mode only - in August 1996 offers the
possibility for a direct comparison of the light-curve with our
simulation. We can reproduce the light-curve, but only for higher
effective temperature than usually assumed. The wavelength dependent
amplitudes (chromatic amplitudes) of this object are not well
reproduced by our simulations, and various possible explanations are
discussed.
Key words: stars: white dwarfs - stars: atmospheres - stars: interiors - convection - stars: variables: general
Variable white dwarfs are non-radial pulsators in which buoyancy
in the gravitational field is the main restoring force of the motions
(so-called g-modes; Chanmugam 1972; Warner &
Robinson 1972) unlike most other variable stars were the
pressure has this role (p- or acoustic modes). The photometric
variability for these stars is a consequence of temperature
fluctuations on the surface of the star (Robinson et al. 1982
= RKN) and not of real geometric changes. The fundamental driving
mechanism is still a matter of some debate.
Brickhill (1991a) concluded that white dwarf pulsations are
driven by a mechanism he called convective driving instead of the
-mechanism usually assumed for most other pulsating
stars. Other authors, however, consider the term "convective
driving'' misleading (see Brickhill 1991a; Kepler et al. 2003; Kotak et al. 2002a and van Kerkwijk et al. 2000 for a more detailed discussion about the driving
mechanism).
The prerequisite for the success of asteroseismology is the correct identification of the radial order (k), the spherical degree (l) and the azimuthal order (m) of the spherical harmonic of each pulsation mode of the star. For stars with many modes present, the standard method is to use the lifting of a degeneracy by rotation or magnetic fields, leading to triplets (for l=1) or quintuplets (l=2) etc., or the use of asymptotic formulae for the period spacing of modes with neighboring k values.
Robinson et al. (1982) developed an alternative method especially useful for pulsators with few modes. This method - further-on called the RKN method - uses the wavelength dependence of the pulsation amplitude, caused by the changing center-to-limb variation. Using a proper normalization, they showed that changes with wavelength should depend only on the "quantum number'' l and not on the overall amplitude or geometric factors as the inclination angle itowards the line-of-sight of the observer. While this method seems ingenuous and straightforward, its application in practice has so far had only limited success.
A central assumption of the RKN method is that the distribution of flux amplitudes on the surface is completely described by spherical harmonics. While this is certainly true for variations in deeper layers of the star, where the pulsations are very nearly adiabatic, it is conceivable that non-adiabatic effects in the outer layers might change this behavior. To test this assumption and to analyze the properties of pulsating white dwarfs in a way different from other studies a numerical simulation code was developed by Ising & Koester (2001, IK), following the pioneering work of Brickhill (1992b and earlier papers) and applied to DAV. They showed that for larger amplitudes and higher spherical degree l the RKN assumption may not always be justified.
In the present study we briefly introduce this numerical model and some minor refinements of it (Sect. 2). In the following sections the study is extended to variable DBs (DBVs) and applied for a comparison of simulations and observations of the prototype of this class, GD 358. We start with the simulation of light-curves of DBVs (Sect. 3), followed by the wavelength-dependent pulsation amplitude spectra (so-called chromatic amplitudes, Sect. 5). Section 4 addresses the question of the determination of the effective temperature in these stars.
The bottom of the model column is at a point sufficiently below the surface convection zone, and deeper than the convection zone may get during a pulsation cycle. The upper boundary is at an optical depth of 10-100 where radiative transfer can be approximated with the diffusion approximation, and within the range of a grid of model atmospheres used to determine the relation between optical depth, temperature and pressure at the outer boundary. Within the column the energy flux has a radiative (described with the diffusion approximation) and a convective component (described with the mixing-length theory, usually with a standard set of parameters like ML2/0.6, which means the set of constants is from the version ML2 and the mixing length is 0.6 in units of pressure scale height; see Tassoul et al. 1990 and Koester et al. 2001 for the nomenclature).
Starting from a static model a periodic sinusoidal pressure and energy
flux variation is superposed at the bottom of the column. As Brickhill
assumed, and Goldreich & Wu (1999a) confirmed, the relative
pressure variation can be taken as constant through the whole
column. The resulting flux variation leads to flux differences
between adjacent cells in the column, which in turn determines the
new temperature stratification
"Observable'' light-curves and chromatic amplitudes are obtained from this local calculation in the following manner. At the lower boundary below the convection zone we assume a distribution of pressure amplitudes a for the whole star according to a spherical harmonic with indices l and m. The whole outer layer is described as a large number of columns, with a variation of the flux at the top according to the specific amplitude at the bottom of that column. The surface flux is then integrated over the visible stellar disk for the observer, taking into account the angle-dependency of the surface intensity ("limb darkening'') and the inclination between the line-of-sight of the observer and the pulsation axis. This results in a flux or spectrum as a function of time, equivalent to an observed light-curve. As a measure of the distortion of the variation (which is purely sinusoidal at the bottom) by the non-linear reaction of the outer layers we use a Fourier decomposition of this light-curve, which in the case of an exactly linear reaction would result in the recovery of one frequency only, but with additional harmonics otherwise. The same analysis is also used for the "local'' calculation, using one column only. The amplitudes of this decomposition are called Fourier coefficients in the following, with the designation "fundamental'' used for the amplitude of the real physical mode, the exciting variation at the lower boundary.
In their work devoted to a study of the principal differences between
this approach and the analytical RKN method, IK used only one stellar
model for a variable DA with
= 11 350 K and
= 8 and a
pulsation period of 100 s. A detailed description of the numerical
model, the exact derivation of Eq. (1), and the major
results can be found in that paper.
The numerical method described in Sect. 2.1 works well
for a large range of parameters as presented in IK. Further studies
showed that very small pressure amplitudes may be affected by
numerical noise, as can be seen easily in Eq. (1). For small
amplitudes the flux differences
become very small compared
to the flux values themselves. As a result the Fourier coefficients
and phases of the emergent flux relative to the input variation show
unphysical erratic behavior.
As a test for the numerical accuracy of the IK-model at small
amplitudes Eq. (1) was linearized with respect to the flux
differences. The fluxes near a static equilibrium solution with
constant flux throughout the column were expanded and first order
linear variations with respect to small temperature changes in the
cell determined. The resulting modified differential equation
for the time-dependent change of the temperature becomes
![]() |
(3) |
![]() |
Figure 1: Relative net flux in the top (#1) and bottom (#31) cell for 0.1% pressure amplitude in a column of for a DAV model. Compared are the results using the standard IK numerical scheme and the linearization described here. Especially at the bottom the numerical noise is reduced significantly (the flux in the bottom cell is multiplied by a factor of 4 for better visibility). |
| Open with DEXTER | |
Higher temperatures and different composition are not the only changes
when simulating DBVs. Whereas the increase in depth of the convection
zone is rather abrupt and steep when going towards lower effective
temperatures in DAV, where a change of 1000 K may switch between off
and strongly on for convection, changes are much smoother in the
DBVs. Even a variation of 4000 K for a 25 000 K model will not switch
off convection completely. IK found a linear regime in the DAV with
sinusoidal light-curves for pressure amplitudes up to about 2%. Above
this value deviations from sinusoidal variations and nonlinear effects
appear, until at about 15% saturation of flux amplitudes occurs and
the light-curves change back to sinusoidal again for higher pressure
amplitudes. In the case of a 25 000 K DBV the sequence of changes is
the same, but the corresponding transitions occur at higher amplitudes
of
10 and
35% pressure amplitude; for a 27 000 K DBV the
values are even higher. This is shown in Fig. 2, which
demonstrates the relation between the pressure amplitude and the
Fourier coefficients of the surface flux for one column in a DBV
model.
Our test cases for simulation of DBVs are models with effective
temperatures of 25 000 and 27 000 K, with a pulsation period of
423 s. The composition is pure He, the input physics, including the
version of mixing length for convection, and the numerical procedures
are identical to those described in IK.
![]() |
Figure 2: Top: Fourier coefficients of the surface flux for one column as a function of pressure amplitude at the bottom. Effective temperature for this model is 25 000 K, the pulsation period is 423 s. Bottom: the same for a 27 000 K DBV model. |
| Open with DEXTER | |
Figure 4 demonstrates the effects of changing the spherical degree l of the mode. From this simulation it is obvious that the light-curve alone does not provide enough information to distinguish between a l=1 and l=2 mode.
![]() |
Figure 3:
Light-curves of a
|
| Open with DEXTER | |
![]() |
Figure 4:
Comparison of a l=1 light-curve (dashed line) with an
inclination of
|
| Open with DEXTER | |
![]() |
Figure 5: Light-curve of GD 358 (full line) observed by Kepler et al. (2003) on August 16th 1996 (upper) and on August 14th 1996 (lower). |
| Open with DEXTER | |
![]() |
Figure 6:
Both: Part of the observed light-curve of GD 358 in August
1996 (dashed line; Kepler et al. 2003). On the upper panel a
fit with a model with 25 000 K equilibrium
|
| Open with DEXTER | |
![]() |
Figure 7: Comparison of August 1996 light-curve (full line) with a sine function (dashed line) of 423 s period. |
| Open with DEXTER | |
![]() |
Figure 8:
Part of the observed light-curve of GD 358 in August 1996
(full line; Kepler et al. 2003) compared with a model with 25 000 K
equilibrium
|
| Open with DEXTER | |
During the last 10 years a large effort was spent in order to solve
the puzzle of the DBV GD 358, the prototype of the class. Three Whole
Earth Telescope (WET) runs (Winget et al. 1994; Vuille et
al. 2000; Kepler et al. 2003) and numerous
single site observations provided a wealth of information on the
pulsation properties. Optical and UV spectra complemented the
photometric observations, but there are still doubts left about such
fundamental parameters as e.g.
.
While almost all optical
observations point toward
=
K, the HST spectra
(Provencal et al. 2000), chromatic amplitudes derived from
Keck observations (Kotak et al. 2002b) and the analysis of
the light-curves in this paper hint at a significantly higher
temperature of about 27 000 K.
In the first WET runs a large number of pulsation modes with periods between 400 and 800 s were identified. Some of the modes found in 1990 could not be re-identified in the 1994 observations, or had noticeable different amplitudes or small frequency shifts, indicating that this object may be much more complicated than expected from theory.
Like all DBVs GD 358 is normally a multi-periodic oscillator and a direct comparison between simulated and observed light-curves is not possible since our numerical code cannot yet calculate several modes simultaneously. However, a very special event for the brightest DBV (V = 13.65) came to our help and provided us with the possibility of a detailed comparison.
In August 1996 Kepler et al. (2003) obtained 48 hours of photometric observations. The light-curve and its Fourier transformation showed a very unusual behavior - only one instead of the usual dozen or so pulsation modes was active and showed a sinusoidal curve with an extreme amplitude of about 44% (peak to peak). Furthermore this single mode was not the usually dominant one with a period of about 770 s but a normally rather weak mode with 423 s. In September of the same year the usual light-curve with about a dozen fundamental modes and a maximum amplitude of 5% was observed again. The cause for this behavior is not understood (see Kepler et al. 2003 for a possible scenario). There is no precedent and it is currently unknown if this is a recurring phenomenon and can happen also in other DBVs. Figure 5 shows the light-curve of GD 358 on two days (16th and 14th) in August 1996.
While the upper panel shows a "normal'' light-curve on the lower panel almost all pulsation energy was dumped into a single mode. The sudden change in the pulsational structure of GD 358 provided the possibility to study this remarkable star during a period of nearly single-mode pulsation, which can be simulated with our numerical procedure.
Figure 6 shows that the simulation is able to reproduce
the light-curve but only using a high
model (but see also
Fig. 8) and a surprisingly high pressure amplitude
a at the bottom of the simulated column. The growing shift between
simulation and observation in both figures and the modulation of the
observed amplitude is probably due to a beating between two very close
excited frequencies. The very nearly sinusoidal form of the light-curve
(see Fig. 7) even with a pressure amplitude of 30% in
the 27 000 K model, and the stronger deviations at 25 000 K are easily
understood from the results presented in Fig. 2 for a
single column.
The improved agreement for higher
does not necessarily mean that
this is the correct temperature, since other changes of the input
parameters can produce similar fits to the light-curve. An example is
given in Fig. 8 where we have used a different
parameterization of the mixing-length theory (ML1/0.5), which
describes convective energy transport with much lower efficiency.
This model gives a fit of similar quality as seen in the lower panel of
Fig. 6 with
= 25 000 K, more in accord with
spectroscopic determinations. The driving amplitudes below the
convection zone are still very high with a=25%.
In the spectroscopic study of Beauchamp et al. (1999) the
authors find better agreement between (time-integrated) observed
spectra and their static models for DBVs (and DABVs) with more efficient convection (ML2/1.25) than we use in our
standard model, which seems to contradict our above finding.
One has to bear in mind, though, that time-dependent observables, even
when averaged over period, may not always be identical to those
quantities in a static model. As an example, for the large amplitude
variations of these models we find that the time-averaged depth of the
convection zone is significantly less than that of a static model with
the same
and the same mixing-length parameters. That brings up
the more general question of the definition and/or measurement of
effective temperatures in stars with flux varying by several 10%.
Figure 9 shows theoretical optical spectra for the
moments of maximum and minimum total fluxes from the integrated disk
for models with equilibrium
of 25 000 and 27 000 K. The spectra
clearly show significant variations over one period. Also shown are
the time-averaged spectra, compared with static models of the same
.
The latter two are almost identical and we conclude that even
at these very large amplitudes the interpretation of (correctly
integrated) spectra with static atmosphere models should result in
correct
estimates. The reason for this is that for the high
temperatures of the DBVs the optical range falls on the Rayleigh-Jeans
part of the source function, where the change of the flux is a linear
function of temperature to first order. This is still true in the UV
range covered by HST spectra; we should expect noticeable
deviations from this only below
1100 Å.
![]() |
Figure 9:
Spectra at maximum (full line) and minimum (long-dashed line)
|
| Open with DEXTER | |
In general these attempts have not been as successful as expected; the agreement of the observed chromatic amplitude with theoretical predictions has often not been convincing. While the numerical IK-model needs its own share of simplifications these numerical integrations avoid central assumptions of the RKN method: that the flux response at the surface is a linear function of the sinusoidal pressure variation below, and that the phase between the two is constant. This is practically equivalent to the assumption, that the spatial distribution of surface temperatures is at all times described by a spherical harmonic function of some degree l.
IK showed that this assumption is indeed often violated, especially at larger amplitude. However, the effects (i and a dependence of chromatic amplitudes) is negligible for l=1, and often small for l=2. Only for higher degrees significant differences from RKN methods should be expected. In this paper we extend that research to the study of fractional amplitudes in DBVs for the special case of GD 358.
![]() |
Figure 10: Observed chromatic amplitude spectra in comparison with a RKN model of l=1 (long-dashed) and l=2 (short-dashed). |
| Open with DEXTER | |
![]() |
Figure 11:
Observed chromatic amplitude spectra in comparison with a numerical
model (long-dashed) and a RKN model (short-dashed). Both models have
|
| Open with DEXTER | |
There are a number of possible explanations for this failure, including a wrong limb darkening in the model atmospheres, or in the assumptions of the RKN methods. Since the model atmospheres describe the time-averaged spectra quite well, as well as the non-pulsating DAs and DBs, which of course depend on the temperature structure as does the limb darkening, we do not expect large errors, unless the variable white dwarfs are somehow fundamentally different from their cousins a few hundred K hotter or cooler. We therefore explore here the possibility that the failure is due to assumptions of the interpretation method, using our completely numerical scheme to explore the parameter space of different amplitudes, effective temperatures, and inclination angles.
![]() |
Figure 12:
Effects of
|
| Open with DEXTER | |
![]() |
Figure 13: The three best fitting models. The observations are shown with a fully drawn line while the fits are the dashed lines. Top: high temperature (28 000 K) fit. Middle: l=2 fit. Bottom: 1st harmonic fit. |
| Open with DEXTER | |
For the l=1 modes the dependence on i and a is negligible, but
is very important. This confirms the result found by Kotak et
al. (2002b) with their RKN-method, who also obtained a much
better fit for higher
(e.g. 27 000 K). With the current
assumptions (
,
,
l = 1, fundamental mode pulsation) for
GD 358 deduced from many different observations we cannot explain the
shape of the chromatic amplitudes.
Apart from using a higher
,
a better fit to the observations is
also possible by assume l=2, or that this mode is a 1st harmonic
(see Fig. 13). In view of the detailed
analysis of Winget et al. (1994), and Vuille et
al. (2000) these possibilities also seem to be unlikely.
Interestingly, however, Kepler et al. (2003) claim the
discovery of a high amplitude l=2 mode with a period of 796 s in the
data from the 2000 WET run. The multiplet structure of this mode is
more complex than the expected quintuplet and may be a hint for mode
interaction. We note that in our study also a l=2 mode (middle panel
in Fig. 13) provides a good fit. Another good fit is
with a first harmonic (lower panel in Fig. 13).
Wu (2001b) predicts larger harmonics for l = 2 modes,
although this may be somewhat uncertain, since already the predictions
for l = 1 harmonics seem to be too large compared to observations.
Nonlinear responses of the surface temperature on pressure variations
in the deeper layers are less important in DBVs than in DAVs, because
the extension of the convection zone does not change as abruptly with
.
The sequence of behavior with increasing pressure amplitude is
similar - from linear to strongly non-sinusoidal and back to simple
sinusoidal light-curves - but the pressure amplitudes for the
transitions need to be significantly higher in DBVs.
The singular event in August 1996, when the prototype star of the DBV
class GD 358 showed a sudden change of the pulsational structure from
a multi-periodic non-sinusoidal small amplitude pulsation to a nearly
single mode with very large amplitude offers a possibility to apply
our simulations directly, which can up to now only simulate one mode
at a time. We are able to simulate the light-curve, and while this
does not help to understand the physical reason for this abnormal
change, it is comforting that we can understand that the light-curve
can still be very nearly sinusoidal at extremely large pressure
amplitudes of 30%. We need a temperature of
= 27 000 K for a
good fit, which is slightly higher than most recent spectroscopic
determinations, although within the range of possible values in the
literature. An alternative solution is a change of the mixing length
parameters, which is equivalent to the assumption of less efficient
convection compared to our standard model. A similar discrepancy is
found for the time-resolved spectroscopy (chromatic amplitudes). The
numerical simulations confirm the simple RKN method, and do not lead
to acceptable fits for the "standard'' values (
,
,
l). Possible solutions for these observations are again a higher
effective temperature (
= 27 000-28 000 K), assuming the mode
studied to be a l = 2 mode or a second harmonic. All these
alternatives lead to significantly better agreement.
One conclusion is that still more work is needed to determine the
basic stellar parameters accurately and to reconcile the different
methods, which result in values between 23 000 and 28 000 K for
.
Is the temperature really at the low end of this range, or
closer to 27 000 K as suggested by our results, as well as some
spectroscopic studies (e.g. Provencal et al. 2000)? The
interpretation with static model atmospheres in the optical and UV
range should lead to correct parameters, but only if proper care is
taken to average correctly over the longest and largest amplitude
periods.
The pulsational properties of this object pose many puzzles. We have mentioned already the August 96 change. A pressure amplitude of 30% (which roughly translates into an optical amplitude of about 60% for the maximum amplitude spot at the surface) is an enormously large variation, far from the regime of linear pulsation analysis at the basis of white dwarf asteroseismology. Is this a very rare event, or can it happen in all pulsators? Has such a sudden mode structure change anything to do with the disappearance of almost all triplets in the 2000 WET data set compared with data obtained in the 1990 and 1994 runs? Also the changes of the pulsation frequencies and their amplitudes between the individual WET runs are in general not understood.
In the 1994 frequency spectrum a large number of third-order harmonics and cross-frequencies were present and affecting some of the major peaks (Vuille et al. 2000), contrary to the results from the 1990 WET run. The spectral observations do not have sufficient temporal resolution to study such detail for 1999 (Kotak et al. 2002b). If such strong mode interactions were present during that observations that might perhaps explain some of the difficulties found in this work and Kotak et al. (2002b). It is unknown at present how the nonlinear interaction between several modes, possibly with additional power near the frequency considered here, would affect the chromatic amplitudes. Also the presence of l=2 modes with similar frequencies may effect the chromatic amplitudes, in particular since for these modes larger amplitudes of the higher harmonics are expected (Wu 2001b). We hope in the future to extend our numerical code to calculate several modes simultaneously present and to study the effect of mode interaction on chromatic amplitudes.
Acknowledgements
This work was supported by the Deutsche Forschungsgemeinschaft, DFG project number Ko 739/18-1. We thank Dr. Rubina Kotak and Dr. Kepler de Souza Oliveira Filho for kindly providing us with their observational results before publication.