Contents

• Numerical simulations of the pulsating DB white dwarf GD 358
• 1 Introduction
• 2 The numerical model
• 3 Light-curve simulations of DBVs
• 4 Temperature estimates from observations
• 5 Chromatic amplitudes
• 6 Conclusions
• References

A&A 406, 657-666 (2003)
DOI: 10.1051/0004-6361:20030739

Numerical simulations of the pulsating DB white dwarf GD 358

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

1 Introduction

Pulsating stars, in particular those with non-radial pulsations and many simultaneously exited modes, provide the possibility to study their interiors by means of asteroseismology, allowing a unique insight into their structure and composition not possible for other stars. In the case of pulsating white dwarfs (DAVs with hydrogen atmospheres and DBVs for helium atmospheres) this not only lead to clues about their structure but in some favorable cases also about their rotation, magnetic fields, interior composition and evolutionary timescales.

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 $\kappa$-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.

   
2 The numerical model

   
2.1 The Ising-Koester model

Based on the work by Brickhill (1992b and earlier) IK developed their extended numerical model. In order to obtain the local time-dependent surface temperature distribution the effect of a flux and pressure variation below the convection zone on the emergent flux is studied. To reduce the numerical effort the model describes not the complete star but a column throughout the outer layers containing the complete convection zone. This simplification is justified as the g-modes (in contrast to the p-modes) propagate predominantly in the envelope of a white dwarf (Hansen & Kawaler 1994).

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 $\Delta F$ between adjacent cells in the column, which in turn determines the new temperature stratification

$$\frac{{\rm d}T}{{\rm d}t} = \frac{\Delta F}{C_{\rm p}~m} + ... ...la_{\rm ad} T \frac{{\rm d}}{{\rm d}t}\frac{\delta P}{P}\cdot$$ (1)

Equation (1) describes the conservation of energy, with $\Delta F$ the difference between the outgoing and incoming flux of a cell. The time-dependent development of the temperature T of a cell of the mass m and the heat capacity $C_{\rm p}$ is calculated by iterating between this equation and the calculation of the radiative and convective fluxes from the resulting new temperatures. At the upper boundary the simulation gives the total flux or effective temperature as a function of time, with the amplitude and period of the pressure variation at the bottom as parameters.

"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 $T_{\rm eff}$ = 11 350 K and $\log g$ = 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.

2.2 Very small amplitudes in the IK-model

The main topic of the present study is the extension of the numerical simulations to variable DB stars. However, before presenting that in the following sections, we turn back to the DAV model used by IK to address some problems noted with the simulations at very small amplitudes. This is an important range, because it should be expected that the results of the numerical analysis should approach the semi-analytical results of Goldreich & Wu (1999a).

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 $\Delta F$ 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

$$\frac{{\rm d}T}{{\rm d}t} = \frac{\Delta F_0 + \delta F}{C_{... ...la_{\rm ad} T \frac{{\rm d}}{{\rm d}t}\frac{\delta P}{P}\cdot$$ (2)

Here $\Delta F_0$ is the flux deviation from an equilibrium model at the beginning of the time step calculated. When constructing a starting model for pulsation calculations, the model is relaxed for a sufficiently long time with constant lower boundary until $\Delta F_0$ is practically zero. $\delta F$ is the small flux change resulting from the time-dependent external pressure excitation $\frac{{\rm d}}{{\rm d}t}\frac{\delta P}{P}$ and calculated by the use of the numerical derivatives

$$\delta F = \left (\frac{\delta F}{\delta T}\right ) \delta T + \left (\frac{\delta F}{\delta P}\right ) \delta P.$$ (3)

Equation (2) is used instead of Eq. (1) for small variations. A comparison of IK results with this linearized version is shown in Fig. 1 for a pressure amplitude of 0.1%. Especially at the bottom of the column the result for the net flux is much smoother than the direct calculation with Eq. (1). A number of numerical experiments with pressure amplitudes from 0.1 to 1.0% and comparisons between linearization and direct solution lead to the following conclusions

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

• For the fundamental mode the Fourier coefficients for amplitude and phase approach closely the numerical results of the IK method at intermediate amplitudes (0.5-0.6%). These coefficients can therefore be extended from the validity range of the IK calculations down to practically zero amplitudes.

• The phase of the fundamental mode is constant and close to the value of $70^{\circ}$ predicted by Wu (1998).

• The variations become very accurately sinusoidal at small amplitudes. Even with the linearization the Fourier coefficients for amplitudes and phases of the higher harmonics cannot be accurately determined; the results depend strongly on minor details of the calculations. It is then preferable to use the predictions of the perturbation analysis if the calculations need to be extrapolated to very small pressure amplitudes.

• The linearization method differs from the IK calculations significantly already at 1% amplitude. The reason is that even at this small amplitude the depth of the convection zone changes considerably in the IK calculation, an effect which is not considered in the linearization.

   
3 Light-curve simulations of DBVs

While the IK paper was focused on DA variables in this paper we extend the study to pulsating helium white dwarfs. Compared to the instability strip for ZZ Ceti between 11 000 and 12 000 K the temperatures for the DBV instability strip are much higher due to the larger temperatures needed to ionize helium. Depending on the model atmospheres used empirically the DBVs are found between 22 000 K and 28 000 K when using pure He atmospheres and between 22 000 K and 25 000 K when traces of hydrogen (not seen directly) are assumed (Beauchamp et al. 1999). For our numerical study we have used pure He atmospheres.

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 $\sim$10 and $\sim$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.

The intensity of the disk was numerically integrated for a spectral range from 1000 to 6600 Å, but for the comparison with observations this was integrated over the Johnson V passband. Calculations were performed until a stationary state was reached. The effect of different pressure amplitudes aexciting the pulsation at the bottom of the considered layer, and of different inclination angles for the observer are shown in Fig. 3. As in the case of the DAV with increasing amplitudes the ascending part of the light curve becomes steeper than the descending part, although this happens only at higher amplitudes as mentioned above.

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 $T_{\rm eff}$ = 27 000 K IK-model. Top: different inclinations i from $i=0^{\circ }$ to $i=90^{\circ }$ in steps of 15 degrees. Bottom: pressure amplitudes a from $a=50\%$ to $a=5\%$ in steps of 10 resp. 5%. The pulsation period is 423 s.

Figure 4: Comparison of a l=1 light-curve (dashed line) with an inclination of $i=60^{\circ }$ with and l=2 at $i=0^{\circ }$ (full line). Both models have an $T_{\rm eff}$ of 25 000 K and a pressure amplitude of $a=40\%$.

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

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 $T_{\rm eff}$ (full line; pressure amplitude $a = 35\%$, inclination $i = 0^{\circ},~l=1$) and on the lower panel with a 27 000 K model (full line; $a=30\%,~i=0^{\circ},~l=1$).

Figure 7: Comparison of August 1996 light-curve (full line) with a sine function (dashed line) of 423 s period.

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 $T_{\rm eff}$ (dashed line; pressure amplitude $a = 25\%$, inclination $i = 0^{\circ},~l=1$) but a different mixing-length theory (ML1/0.5) than used for the upper panel of Fig. 6.

   
3.1 GD 358 in August 1996

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.  $T_{\rm eff}$. While almost all optical observations point toward $T_{\rm eff}$ =  $24~000\pm 1000$ 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 $T_{\rm eff}$ 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 $T_{\rm eff}$ 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 $T_{\rm eff}$ = 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 $T_{\rm eff}$ 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%.

   
4 Temperature estimates from observations

For an assumed pressure amplitude of a=30%, $T_{\rm eff}$ as measured through the total energy output at the surface of one local column varies between 22 700 and 31 400 K for the model with equilibrium $T_{\rm eff}$ of 27 000 K. These are also the extremes of the variation on the surface. The corresponding numbers for the 25 000 K equilibrium model are 21 800 and 28 700 K. In this section we study the effect of this enormous variation on the empirical determination of $T_{\rm eff}$ from time-averaged spectra.

Figure 9 shows theoretical optical spectra for the moments of maximum and minimum total fluxes from the integrated disk for models with equilibrium $T_{\rm eff}$ 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 $T_{\rm eff}$. 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 $T_{\rm eff}$ 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 $\sim$1100 Å.

Figure 9: Spectra at maximum (full line) and minimum (long-dashed line) $T_{\rm eff}$ for a model with 25 000 K equilibrium $T_{\rm eff}$ (upper) and 27 000 K (lower), compared to non-pulsating models (short-dashed line) and time-averaged (dotted line) of the same $T_{\rm eff}$.

   
5 Chromatic amplitudes

The wavelength dependence of pulsation amplitudes has in recent years been established as an alternative method to identify pulsation modes in variable white dwarfs, in particular in the ZZ Ceti with only few modes excited (Robinson et al. 1995; Clemens et al. 2000; Kotak et al. 2002a). It has even been possible to identify surface motions, believed to be associated with the pulsations, using time resolved spectra (van Kerkwijk et al. 2000). Robinson et al. (1982) have laid the foundation for all these studies by demonstrating that after proper normalization of the amplitudes they become independent of the inclination of the pulsation axis to the line of sight i and the amplitude of the pulsation a, and depend solely on the spherical degree of the mode l. The physical basis of this effect is the wavelength-dependent limb-darkening of the stellar atmosphere: the cancellation effect of opposite variations on the surface becomes less important the stronger the limb-darkening, e.g. towards the UV end of the spectrum. Almost all recent attempts to use this method follow the work of RKN with some minor improvements and thus share the approximations and implicit assumptions of that result.

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.

5.1 Pulsation modes of GD 358

The GD 358 WET runs lead to the identification of more than hundred pulsation modes in this star. A detailed analysis of these modes made it further possible to identify them with l=1 modes of different radial node number, which have different amplitude distributions with depth and therefore can - at least in principle - sample different regions inside the star. The triplet structure of the modes clearly points toward l=1 as expected for a rotational splitting in 2l+1components, and the splitting of 39.2 s is in agreement with the theoretical predictions (Winget et al. 1994).

Figure 10: Observed chromatic amplitude spectra in comparison with a RKN model of l=1 (long-dashed) and l=2 (short-dashed).

In June 1999 GD 358 was observed by Kotak et al. (2002b) with high resolution time-resolved spectroscopy at the Keck II telescope. The observations and analysis concerning identified modes and velocity fields are described in Kotak et al. (2002b). Figure 10 shows the chromatic amplitudes for the dominant mode at 776.42 s. They are compared with theoretical predictions using our own version of an RKN method for modes with l=1 and l=2. Here $T_{\rm eff}$ = 24 000 K is assumed - following the analysis of WET runs for this star. The model does not predict correctly the shape of the amplitudes in the broad wings of the HeI lines, and the overall agreement is quite unsatisfactory.

Figure 11: Observed chromatic amplitude spectra in comparison with a numerical model (long-dashed) and a RKN model (short-dashed). Both models have $T_{\rm eff}$ = 24 000 K, $\log g$ = 8.0 and l=1. The pressure amplitude for the numerical simulation is $a=20\%$ and the inclination $i=60^{\circ }$.

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 $T_{\rm eff}$ , pressure amplitude a and spherical degree l on the chromatic amplitude spectra. The uppermost panel shows a $T_{\rm eff}$ sequence for l=1 and $a=20\%$ the second the same for l=2 and the lower two panels for l=1, 2 and $a=5\%$.

Figure 11 shows the direct comparison of both methods using $T_{\rm eff}$ and mode number l as obtained from observations. It is obvious that both theoretical methods give nearly identical results, but that neither of them fits the observations. The curvature between the central cores of adjacent HeI lines in particular is not correctly represented by the models. In Fig. 12 we explore the effect of changing the parameters (a, l, $T_{\rm eff}$). The results of this parameter study can be summarized as follows.

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.

For the l=1 modes the dependence on i and a is negligible, but $T_{\rm eff}$ 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 $T_{\rm eff}$ (e.g. 27 000 K). With the current assumptions ( $T_{\rm eff}$, $\log g$, 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 $T_{\rm eff}$, 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.

   
6 Conclusions

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 $T_{\rm eff}$. 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 $T_{\rm eff}$ = 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 ( $T_{\rm eff}$, $\log g$, l). Possible solutions for these observations are again a higher effective temperature ( $T_{\rm eff}$ = 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 $T_{\rm eff}$. 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.

References

 
• Beauchamp, A., Wesemael, F., Bergeron, P., et al. 1999, ApJ, 516, 887 NASA ADS  
• Bradley, P. A., & Winget, D. E. 1994, ApJ, 430, 850 NASA ADS  
• Bradley, P. A. 2001, ApJ, 552, 326 NASA ADS  
• Brickhill, A. J. 1983, MNRAS, 204, 537 NASA ADS  
• Brickhill, A. J. 1990, MNRAS, 246, 510 NASA ADS  
• Brickhill, A. J. 1991, MNRAS, 251, 673 NASA ADS  
• Brickhill, A. J. 1991, MNRAS, 252, 334 NASA ADS  
• Brickhill, A. J. 1992, MNRAS, 259, 519 NASA ADS  
• Brickhill, A. J. 1992, MNRAS, 259, 529 NASA ADS  
• Chanmugam, G. 1972, Nat. Phys. Sci., 236, 83 NASA ADS  
• Clemens, J. C., van Kerkwijk, M. H., & Wu, Y. 2000, MNRAS, 314, 220 NASA ADS  
• Goldreich, P., & Wu, Y. 1999a, ApJ, 511, 904 NASA ADS  
• Goldreich, P., & Wu, Y. 1999b, ApJ, 523, 805 NASA ADS  
• Hansen, C. J., & Kawaler, S. D. 1994, Stellar Interiors (New York, Inc.: Springer-Verlag)  
• Ising, J., & Koester, D. 2001, A&A, 374, 116 NASA ADS  
• Kepler, S. O., Robinson, E. L., Koester, D., et al. 2000, ApJ, 539, 379 NASA ADS  
• Kepler, S. O., Nather, R. E., Winget, D. E., et al. 2003, A&A, 401, 639 NASA ADS  
• Koester, D., Allard, N. F., & Vauclair, G. 2001, A&A, 291, L9 NASA ADS  
• Kotak, R., van Kerkwijk, M. H., & Clemens, J. C. 2002, A&A, 388, 219 NASA ADS  
• Kotak, R., van Kerkwijk, M. H., Clemens, J. C., & Koester, D. 2003, A&A, 397, 1043 NASA ADS  
• Provencal, J. L., Shipman, H. L., Thejll, P., et al. 2000, ApJ, 542, 1041 NASA ADS  
• Robinson, E. L., Kepler, S. O., & Nather, R. E. 1982, ApJ, 259, 219 NASA ADS  
• Robinson, E. L., Mailloux, T. M., Zhang, E., et al. 1995, ApJ, 438, 908 NASA ADS  
• Tassoul, M., Fontaine, G., & Winget, D. E. 1990, ApJS, 72, 335 NASA ADS  
• van Kerkwijk, M. H., Clemens, J. C., & Wu, Y. 2000, MNRAS, 314, 209 NASA ADS  
• Vuille, F., O'Donoghue, D., Buckley, D. A. H., et al. 2000, MNRAS, 314, 689 NASA ADS  
• Warner, B., & Robinson, E. L. 1972, Nat. Phys. Sci., 234, 2  
• Winget, D. E., Nather, R. E., Clemens, J. C., et al. 1994, ApJ, 430, 839 NASA ADS  
• Wu, Y. 1998, Excitation and Saturation of White Dwarf Pulsations, PhD Thesis, California Institute of Technology, Pasadena, California  
• Wu, Y., & Goldreich, P. 1999, ApJ, 519, 783 NASA ADS  
• Wu, Y., & Goldreich, P. 2001a, ApJ, 546, 469 NASA ADS  
• Wu, Y. 2001b, MNRAS, 323, 248 NASA ADS